mMDS. The Shepard diagram for the present case ('Graph3'), demonstrating both linearity and the perfect suitability of a zero intercept, is shown below.*
[![10._Shepard_plot_mussels_mMDS_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/10-shepard-plot-mussels-mmds-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/10-shepard-plot-mussels-mmds-i.png)
# 13.4 Example: Gulf of Maine invertebrates - functional resemblance
There are many situations where the standardisation of samples is required as a pre-treatment prior to analysis, but which needs to be done ***separately within groups*** of variables that may be identified by an indicator.¶ Here, we shall consider a dataset comprised of occurrences of $p$ = 91 macroinvertebrate species recorded from intertidal areas at each of $N$ = 12 exposed rocky headlands (surveyed in 2012) in the Gulf of Maine (Fig. 13.3, {{@954#bkmrk-trott2022}}†). For each species, we have additional information in the form of trait data. The trait data for each species consists of presence/absence information for a series of $q$ = 93 individual trait variables. The trait variables are organised into 14 groups. These groups of traits are: Form, Position, Lifestyle, Trophic Relation, Body Shape, Body Support, Flexibility, Ecoengineer, Fertilization, Development, Size, Body Plan, Asexual Reproduction and Regeneration. The number of categories (trait variables) within each trait group varies, and any given species can be recorded as a '1' (indicating that they belong) to one or more of these categories within each trait group.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/15b-gulf-of-maine-google-map-mark-up.png)
*Fig. 13.3 Map of the Gulf of Maine showing 12 sites where macroinvertebrates were sampled by {{@954#bkmrk-trott2022}}. Satellite image: Google Earth.*
Our initial agenda here may be to consider the relationships among the sites based on the species they contain in the usual way, e.g., using the Bray-Curtis (or Sørensen) resemblance measure directly on the species $\times$ site matrix. However, we may wish to nuance this calculation further by incorporating relationships among the species, based on their traits. In other words, if two sites do not share any species in common, they might nevertheless share two species that have similar traits. We can exploit the method of calculating taxonomic resemblance ('**Gamma+**', {{@954#bkmrk-clarkeetal2006b}}) in order to calculate, instead, a ***functional resemblance***, using a ***matrix of inter-species relationships built from the trait data*** (e.g., {{@954#bkmrk-myersetal2021}}).
Further to this aim, we shall consider our trait data matrix such that the classical role of 'samples' is given to the species (among which we shall calculate the resemblances), while traits are given the role of 'variables'. However, we want each trait ***group*** to be weighted equally in the calculation. To make sure the sum of values for each species within each trait group sums to 100 (so gets equal weight), we will want to apply a pre-treatment to standardise across the trait variables for each species ('sample'), but separately within each trait grouping. This standardisation would not be necessary if every species only had one trait within a group, but that is not the case here. For example, the mollusc species, *Adalaria proxima* (abbreviated ADP) is in two trophic categories: it grazes on algal fronds and blades (Tr-P-GF), and also grazes macro prey on the substratum (Tr-P-GSM).
Our analysis pathway looks like this:
- Open the trait data matrix in PRIMER 8.
- Standardise the 'samples' (species here) across the trait variables, separately within each trait grouping, using our new tool in P8, **Pre-treatment** > **Standardise...**.
- Calculate functional resemblances among species using Bray-Curtis, based on the standardised trait data.
- Examine the inter-species functional relationships in an ordination (nMDS).
- Calculate functional resemblances among the 12 sites, by incorporating these species inter-relationships.
- Examine the functional relationships among the sites in an ordination (nMDS).
#### Open the trait data matrix
1. In PRIMER 8, click **File** > **Open...** and open up the file named 'Gulf_of_Maine_invert_traits.pri' (found in the 'Gulf_of_Maine_inverts' folder inside the 'Examples_P8' folder), as shown below.
[![16._GoM_Traits_data_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/16-gom-traits-data-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/16-gom-traits-data-i.png)
Note here that the 'Samples' are actually the species (the names are abbreviated). More information about the species can be seen by clicking **Edit** > **Factors...**. Note also that the traits are variables. More information about the traits (also abbreviated) can be seen by clicking **Edit** > **Indicators...**.
The first trait group is actually the name of the Phylum for each species. These are not actually traits, so let's start by selecting all of the traits that do not belong to the 'Phylum' group.
2. From the 'Gulf_of_Maine_invert_traits' matrix, click **Select** > **Variables...**, then choose
($\bullet$Indicator levels Trait category) and click the 'Levels...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/17e-levels-button.png)). In the Selection dialog, first click the double right arrows button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/17d-double-right-arrow.png)) to move all of the groups to the 'Include' column on the right. Then click on 'Phylum' and the left arrow button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/17c-left-arrow.png)), so that it appears back in the 'Available' column on the left (as shown below), then click '**OK**' (in both dialog windows).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/17-select-variables-gom.png)
This will create a dataset that is re-coloured with a blue background, showing only the selected subset of traits (excluding 'Phylum'), like this:
[![18._Subsetted_traits_GoM_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/18-subsetted-traits-gom-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/18-subsetted-traits-gom-i.png)
Any analyses from this matrix will now only be done on the sub-setted data. To keep things tidy, click **Tools** > **Duplicate** and the subsetted data will now be provided in the Explorer tree as a new matrix called 'Data1'. We can re-name this to 'Traits' by clicking **File** > **Rename Data** and typing in this desired new name, then clicking '**OK**', as shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/19-rename-dialog.png)
#### Standardise the trait data
3. From the 'Traits' data sheet, click **Pre-treatment** > **Standardise...** and choose (Standardise $\bullet$Samples) & (By $\bullet$Total) & (Standardise within groups/levels of $\checkmark$Indicator Trait category) & (Output $\bullet$Percentages), then click '**OK**', as shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/20-standardise-within-groups-dialog.png)
The resulting data sheet (called 'Data1') can be renamed (using **File** > **Rename Data**, as we have done before) to (say) 'Standardised Data', for clarity. The standardised trait data sheet will look like this:
[![21._Standardised_trait_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/21-standardised-trait-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/21-standardised-trait-i.png)
Note how the distribution of traits within a species for any particular trait group (in this example, the trait groups are identifiable by reference to the letters appearing before the dash '-' in the variable names, e.g., 'Po', 'Li', 'Tr', etc.) will sum to 100.
#### Calculate relationships among species based on traits
4. From the 'Standardised Data' sheet, click **Analyse** > **Resemblance...** > (Measure $\bullet$Bray-Curtis similarity) & (Analyse between $\bullet$Samples), as shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/22-resemblance-dialog.png)
The resulting resemblances among the species, based on their traits is shown below.
[![23._Resem_among_species_'samples'_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/23-resem-among-species-samples-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/23-resem-among-species-samples-i.png)
We mustn't forget that the role of 'samples' is being taken by the species here! (The traits were the variables that created this matrix). However, in subsequent analyses to follow, the species will be variables, and sites will be samples. So let's go ahead and change the properties of this resemblance matrix now.
From 'Resem1', click **Edit** > **Properties...** and change the following:
- Title 'Gulf of Maine invertebrates'; and
- Between $\bullet$Variables
then click '**OK**' (as shown below).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/23b-resem-among-species-variables.png)
We now have a resemblance matrix among the species which are clearly identified as variables, hence that is ready for ensuing analyses (see below).
[![23c._Resem_among_species_'variables'_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/23c-resem-among-species-variables-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/23c-resem-among-species-variables-i.png)
As a next step, let's aim to visualise similarities among species, based on their traits, using a non-metric MDS ordination.
#### Visualise inter-species trait-based relationships
5. From 'Resem1', click **Analyse** > **MDS** > **Non-metric MDS (nMDS)...**, and go with all of the default options in the dialog (just click '**OK**'). The resulting 2D ordination ('Graph1' in the Explorer tree) is shown below. Different symbols/colours show the different phyla in which individual species (labeled by their abbreviations) belong.
[![24._Trait-space_nMDS_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/24-trait-space-nmds-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/24-trait-space-nmds-i.png)
This ordination has rather high stress (0.214), so we might opt to examine the 3D graphic ('Graph1' in the Explorer tree), which has a lower stress (0.138), as shown below. In essence, this ordination shows us the relationships among these invertebrate species in multi-dimensional trait space.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/25-3d-nmds-traits.gif)
Next, we shall use these trait-based relationships among the species to inform the analysis of functional turnover among the sites, based on the species they contain.
#### Calculate functional resemblances (Gamma+) among sites
6. First, we need to open up the species $\times$ site matrix of presence/absence (occurrence) data in the same workspace (click **File** > **Open...**). This file is called 'Gulf_of_Maine_invert_occurrences.pri' (also found in the 'Gulf_of_Maine_inverts' folder inside the 'Examples_P8' folder).
[![25._Species_by_Site_GoM_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/25-species-by-site-gom-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/25-species-by-site-gom-i.png)
The sites have abbreviated names. Full names (as shown in the map in Fig. 13.3 above) can be seen by clicking on **Edit** > **Factors...**. You will also see that the sites have been given a rank ordered value for their position along the Maine coastline ('West to East'), and have also been classified as occurring either north or south of Penobscot Bay, viz:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/26-factors-gom.png)
7. Now we shall calculate functional (trait-based) resemblances (Gamma+) among the sites. From the 'Gulf_of_Maine_invert_occurrences' worksheet, click **Analyse** > **Resemblance...**, and in the 'Resemblance' dialog window:
- Choose (Analyse between $\bullet$Samples) & (Measure $\bullet$Other > $\checkmark$ Taxonomic P/A > Gamma+), and click the 'Taxonomy...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/27b-taxonomy-button.png)).
- In the 'Variable Relationship' dialog window, choose (Type $\bullet$Resemblance), and click the 'Details...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/27d-details-button.png)).
- In the 'Resemblance (Data)' dialog, choose (Variable resemblance worksheet: Resem1).
- Click '**OK**' (3 times, once for each successive dialog window) to complete the operation.
The dialog windows involved in this operation (i.e., to calculate Gamma+ on the basis of a resemblance matrix, which in our case is a trait-based resemblance matrix) are shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/27-functional-resemblance-operation.png)
The resulting resemblance matrix of Gamma+ values among the sites (called 'Resem2') is shown below.
[![28._Gamma+_trait_resem_among_sites_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/28-gamma-trait-resem-among-sites-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/28-gamma-trait-resem-among-sites-i.png)
#### Visualise functional turnover among sites
8. We can now easily do an ordination (or indeed a cluster analysis or any other resemblance-based analysis) of the trait-based turnover among the sites. From 'Resem2', click **Analyse** > **MDS** > **Non-metric MDS (nMDS)**, just keep all of the defaults and click '**OK**'.
9. From the resulting 2D nMDS ordination graphic (called 'Graph5' in the Explorer tree), click **Graph** > **Sample Labels & Symbols...** and choose to display labels by the factor of 'Site name', and symbols by the factor of 'Penobscot Bay', as shown in the dialog below, then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/29-labels-symbols-dialog-gom.png)
The resulting nMDS graphic looks like this:
[![30._nMDS_functional_turnover_nMDS_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/30-nmds-functional-turnover-nmds-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/30-nmds-functional-turnover-nmds-i.png)
This plot shows there is a clear shift in functional traits of intertidal assemblages from sites located to the south (blue symbols) *vs* those located to the north (amber symbols) of Penobscot Bay. We can do statistical tests of relevant hypotheses about this using ANOSIM (from 'Resem2', click **Analyse** > **ANOSIM...**). Doing this shows that:
- This shift is statistically significant (one-way unordered ANOSIM of the factor 'Penobscot Bay', $R$ = 0.546, $P$ = 0.0022); and
- There is a significant sequential change in trait-based relationships among these communities along the coastline (one-way ordered ANOSIM of the factor 'West to East', $R^O$ = 0.422, $P$ = 0.0026).
---
¶*Equally, and in a directly analogous fashion, there are situations where the standardisation of **variables** is required as a pre-treatment prior to analysis, but which needs to be done separately within groups of **samples** that may be identified by a **factor**.*
---
†*See also the following abstract by Thomas J. Trott, highlighted in the 4th World Conference on Marine Biodiversity: ['Traits matter: when rarity means more than abundance to functional diversity'](https://peerj.com/preprints/26840v1/).*
# 14. Create ordered groups
# 14.1 Overview
There is a new tool in PRIMER 8 that permits the end-user to ***create a set of ordered groups***, based on the numerical values of a given variable.
#### Rationale ####
Below we provide a few examples of this tool's potential utility. There are many more!
### Binning and consolidating information ###
A case where we might want 'categories', but where there may be (initially) continuous values, occurs when we have a broad-scale dataset with a *lot of replicates* that span a very *large area* or domain (in terms of latitude and longitude). For example, consider the [Continuous Plankton Recorder (CPR Survey](https://www.cprsurvey.org/). Many other ocenographic, landscape-level and environmental datasets have these characteristics. We might want to ***consolidate*** that information in some way, perhaps in order to make comparisons with more discrete (less continuous) datasets or variables (e.g., dis-continuous biotic samples at particular sites or regions). We could consider making a 'grid' over the sampling extent (i.e., by creating 1-degree sized bins for each of latitude and longitude) and then get the average values of the variable(s) within those grid cells. Indeed, to understand spatial relationships and patterns, it may be essential to bin samples into latitude (and/or longitude) groupings, so that structural changes through space can be more easily summarised, tested, visualised in plots, etc.
### Classification based on an ordination axis ###
Let's consider another scenario. Perhaps you have just done a principal component analysis (PCA) on a set of morphometric variables (e.g., for a set of individuals belonging to a given species). Suppose the first PC axis explains a lot of the multivariate variation and the position of any individual on that PC axis can essentially be interpreted as a measure of its overall relative size. The values of the individuals on the PC axis are continuous, but maybe you want to classify the individuals into (say) five roughly equal-sized groups along that axis (from small to large). In the past, this task could have taken some considerable time, sorting and data wrangling to get your groups based on the PC axis. Our new tool in PRIMER 8 will create equal sample-sized groups along your chosen continuous variable (e.g., a PC axis, or any other meaningful ordination axis of choice) in a flash.
### Setting known break points in the data
Yet another use-case for this tool is the situation where you want to identify a set of samples that fall within a certain range of values for a continuous variable; for example, suppose you want to separate out samples occuring at sites that are below a certain temperature (e.g., 20°C). You can use the tool to specify the value of 20 as a 'break point', and get a factor that groups samples either side of that break point. Multiple such breaks points can also, optionally, be constructed all at once. For example, suppose you have the continuous variable of depth (in metres); you might want to 'bin' samples into several categories according to their depths, using break points of 10, 20, 30, 40, etc. (We shall see an example of this in the next section.)
### Finding natural peaks or gaps to define groupings
You may have a continuous variable that is not evenly distributed across its full range. For example, a histogram may show a multimodal distribution pattern. You may want to define groups that coincide in some natural way with that uneven distribution, such that each group will capture an individual 'mode' or 'peak'. You may not be certain about where the break-points in-between those modes should be placed. In this case, you would want the tool to find suitable break points for you by (say) minimising the within-group sum of squares for the variable, given the number of groups (peaks). If the modes are asymmetric, then sums of distances-to-median (rather than the within-group sum of squares) might be more appropriate to use here. Or you could go even more non-parametric in your criterion and choose groups so as to minimise only the rank inter-point Euclidean distances within each group.
### Specifying groups of replicates for dispersion weighting
Many species abundance variables (e.g., counts) display intrinsic mean-variance relationships. In other words, the variance increases with the mean. This relationship is particularly strong and important to consider in the context of species that tend to aggregate, e.g., such as fish that school, birds that flock, mammals that travel in herds, insects that cluster together in swarms, etc. We may wish to accommodate this by applying a pre-treatment to our data called ***dispersion weighting*** (see {{@954#bkmrk-clarkeetal2006a}} for details). This pre-treatment option effectively transforms our data to consider the abundances (counts) of ***clusters*** of organisms (for those species that do indeed demonstrate significant clustering of individuals), rather than analysing the raw abundance values. Of course, the dispersion weights that may be needed here will vary according to the individual species' distribution of counts across the sampling units.
To characterise mean-variance relationships, hence to apply dispersion-weighting, we need ***groups*** (levels of some factor) across which we can calculate and measure means and variances. In cases where we do not have any groups to begin with (e.g., whenever our sampling units occur along a continuous gradient, in a simple spatial array or regular samples through time, etc.), then there will be no way to use dispersion weighting at all. Usnig this new tool, we can create groups along one (or more) continuous variables which can then form the basis of a dispersion-weighting pre-treatment option.
### Plotting symbols and colours ###
It can be extremely useful to create groups from a continuous variable for rather simple practical or logistic reasons. For example, suppose we want to use different colours or symbols for different levels of a continuous variable, such as temperature. We may not want to have a different colour for every single temperature, but rather we might prefer to have a colour (and/or symbol) for all cases where temperature values occur within some specified range, for each of a set of specified temperature groups (eg., low, medium and high). The use of well-defined symbols and colours like this can make interpretations of patterns in ordinations, especially along gradients, much easier to visualise and communicate to others.
#### Implementation
This tool is accessed by clicking **Tools** > **Create Ordered Groups...**. The dialog window is shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/01-ordered-groups-dialog.png)
There is clearly a large diversity of options to choose from to create ordered groups based on a chosen continuous variable, as shown in the dialog above. The various methods that can be used to create groups are itemised and outlined, correspondingly, in more detail below.
#### Methods for creating groups
Consider an individual variable $Y$, which has a set of observed sample values $y$ and which our algorithm is going to split into $g$ ordered groups using some criterion. Let the number of sample values in each group be denoted by $n_i$ for $i = 1, \ldots, g$ and the total number of sample values is $\sum_i^g n_i = N$. Also, let $y_{ij}$ denote the $j$th observed value of $Y$ in the $i$th group for $j = 1, \ldots, n_i$ and $i = 1, \ldots, g$. In what follows, we shall let $\bar{y}_ i = \sum_{j=1}^{n_i} {y_{ij}} / {n_i} $ be the sample average of the variable for group $i$, and $m_i$ be the median value for the set of samples in group $i$.
The methods for choosing ordered groups based on the variable $Y$ include:
- (1). Specify the ***number of groups***, $g$, and, given that number, generate groups:
- (a). with ***equal sample sizes*** ($n_1 = n_2 = \cdots = n_g$); or
- (b). at ***equally spaced intervals***, so as to ***minimise***:
- (i). ***the within-group sum of squares***, i.e., $\sum_{i=1}^{g}\sum_{j=1}^{n_i} (y_{ij} - \bar{y}_ i)^2$; or
- (ii). ***the sum of within-group mean absolute deviations (MAD) from the median***, i.e., $\sum_{i=1}^{g}\sum_{j=1}^{n_i} |y_{ij} - m_i|$; or
- (iii). ***the sum of within-group average rank inter-point Euclidean distances***.
- (2). Specify ***break(s)*** manually, i.e., a list of specific values of the variable that will serve as break-points between consecutive ordered groups (e.g., 10, 20, 50, 100); or
- (3). Specify the ***quantile(s)*** of the variable's empirical distribution where you want break-points to be (e.g., 0.25, 0.5, 0.75); or
- (4). Specify groups by nominating a certain value for the ***number of samples in each group***, i.e. $n_i = n$ for all $i$.¶
---
†*Dispersion weighting ({{@954#bkmrk-clarkeetal2006a}}) is achieved easily and directly in PRIMER 8 by clicking **Pre-treatment** > **Dispersion Weighting...**.*
---
¶*Note that if there are any 'remainders', given the total sample size, these will be placed in a final group, corresponding to the largest values of the variable, but with reduced $n$. For example, if you have 10 samples and you ask for groups having a sample size of $n$ = 3, the tool will create 3 ordered groups {1,2,3}, {4,5,6}, {7,8,9}, and one remainder sample will spill into a fourth group on its own {10}.*
# 14.2 Example: NE Pacific groundfish vs depth
To demonstrate the creation and use of ordered groups from a continuous variable, we will look at data comprised of an excerpt from the West Coast Groundfish Bottom Trawl (Slope and Shelf Combination) Survey, conducted annually by the National Oceanic and Atmospheric Association (NOAA)'s Northwest Fisheries Science Center¶ and available online ([https://www.nwfsc.noaa.gov/data/map](https://www.nwfsc.noaa.gov/data/map)). This specific excerpt was created and used by {{@954#bkmrk-andersonetal2022}} and is also available on [Dryad](https://datadryad.org/dataset/doi:10.5061/dryad.c59zw3rbp). Data are provided from trawl surveys for the years from 1999-2018 inclusive. There are two data sheets referred to here (both are found in the Examples_P8 > NE_Pacific_groundfish folder), namely: (1) counts of abundances of $p$ = 310 species of groundfish obtained in each trawl (in the file called 'NE_Pacific_groundfish.pri') and (2) values for latitude, longitude, depth (in m) and the area swept per trawl (in the file called 'NE_Pacific_env.pri').
In this example, there are $N$ = 6002 (!) rows of data (sample trawls), with depth values from 24 m - 1,428 m, and latitudes from 32°N to 48°N. Constructing a nMDS ordination plot of these fish assemblages for all 6002 samples is unlikely to be helpful here! To make some sense of these data that span such large latitudinal and depth ranges, it would be helpful to create some groups of samples occuring at similar depths (and also groups of samples occurring at similar latitudes, e.g., see {{@954#bkmrk-andersonetal2013}}). Here, we shall focus purely on creating ordered groups of samples from the variable of ***depth***. More specifically, our interest lies in uncovering potential changes in the structure of fish assemblages along the depth gradient.
#### Open file and create ordered depth groupings
1. Open up the file called 'NE_Pacific_env.pri' in PRIMER 8. It will look like this:
[![02._Groundfish_Env_data_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/02-groundfish-env-data-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/02-groundfish-env-data-i.png)
2. From the NE_Pacific_env data sheet, click **Tools** > **Create Ordered Groups...**, like so:
[![02b._Groundfish_Env_data_click_Tools_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/02b-groundfish-env-data-click-tools-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/02b-groundfish-env-data-click-tools-i.png)
3. In the resulting dialog, choose:
- Variable: Depth (m)
- Grouping criterion > $\bullet$Specify break(s) 50, 100, 200, 400, 600, 800, 1000, 1200
- Output >
- Factor name: Depth bin output
- Factor level labels: $\bullet$Lower bound (LB)
- Output as $\checkmark$Factor
- Output group information to ($\checkmark$Worksheet) & ($\checkmark$Histogram with group boundaries)
The dialog with these choices will look like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/03-dialog-create-groups.png)
#### Output from the 'Create Ordered Groups...' tool
The output file ([![04a._Notepad_icon_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/04a-notepad-icon-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/04a-notepad-icon-i.png)) called 'Ordered Groups1' in the Explorer tree, has some useful information about how the groups were created. For each group we can see: the minimum, median, maximum, lower bound (LB) upper bound (UB), the mean of (LB) and (UB), and the number of samples that fell into each group. This output file is shown below:
[![04._Output_rtf_create_groups_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/04-output-rtf-create-groups-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/04-output-rtf-create-groups-i.png)
(As an aside, note that we opted in this example also to output this information to a separate data sheet as well, which is called 'Data1' in the Explorer tree.)
Another important part of our output (called 'Graph1' in the Explorer tree) is a ***histogram*** of the 'Depth (m)' variable, with the values for the breaks separating the groups shown as vertical dotted lines (as shown below). This is a really helpful way to see the groupings we obtained by reference to the full distribution of sample values for the variable of interest.
[![05._Depth_groupings_hist_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/05-depth-groupings-hist-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/05-depth-groupings-hist-i.png)
Another thing we get (and perhaps the most useful 'handle' for subsequent analyses we may want to do) is a new ***factor*** that is now associated with our original data file. To see this factor, click on the NE_Pacific_env data sheet in the Explorer tree, then click **Edit** > **Factors...**. In the 'Factors' window, you will now see this new factor, called 'Depth bin output', provided as the last column (on the right), like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/06-factors-sheet-a.png)
#### Tweak the factor level names (optional)
We asked for the lower bound value of each group to be made the factor level names of our groups, but we can see that these are not necessarily whole numbers. It might be nice to round these to the appropriate whole number, in each case.
4. First, let's duplicate the newly created factor. From the NE_Pacific_env data sheet, click **Edit** > **Factors...**. In the 'Factors' window, click anywhere in the column named 'Depth bin output', then click the 'Duplicate' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07a-duplicate-button.png)).
5. You will get a new column with a duplicate factor called 'Depth bin output1'. Click anywhere in this new column, then click the 'Rename...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07b-rename-button.png)). Rename this factor simply as 'Depth', then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07c-rename-factor-dialog.png)
6. Now we are going to rename the levels of this factor called 'Depth' to whole numbers. Click the 'Rename Levels...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07d-rename-levels-button.png)).
In the 'Rename Factor Levels' dialog (new to PRIMER8!), you will see two columns: the 'Existing Level Name' on the left and the 'New Level Name' on the right. Change the values in the 'New Level Name' to the desired values , and click'**OK**', as shown below.
[![07e._Rename_levels_sheet_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/07e-rename-levels-sheet-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/07e-rename-levels-sheet-i.png)
In our case, we have specified new factor level names that are still ***numeric*** and consist of whole numbers that will make sense for us in this example for plots/symbols, etc. But note that this tool can be used to re-define the factor level names to anything we wish (not necessarily numbers). Note also, for this example, that we have to be careful not to do too much 'rounding'. We need to stay true to what we know about the data and the bounds of the groups we have created. Bear in mind that we could have changed the names of these factor levels to ranges of depths (e.g., such as 50-100m), which may be better (or more accurate). However, if we change these names to ranges in this way, then we no longer have a strictly numeric factor. Factor levels that are numeric can be really useful in PRIMER, because they allow us to do things like treat the factor as ***ordered*** in an ANOSIM, or superimpose ***trajectories*** to connect consecutive depths on an ordination plot, etc. In this example, given the new names we have chosen, whenever we describe this factor we will have to be clear what the labels mean; specifically, that the group that we have named '50' here corresponds to samples that occurred between 50 m and 100 m in depth, and so on.
#### Create a factor for Latitude
We can repeat the above steps for another important spatial factor here: namely, latitude.
7. From the 'NE_Pacific_env' data sheet in the Explorer tree, click **Tools** > **Create Ordered Groups...**, then choose to create ordered groups from the variable of 'Latitude (dd)', and specify the breaks to occur in 2-degree increments: {34, 36, 38, 40, 42, 44, 46, 48}, as shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/11-latitude-dialog-to-create.png)
A histogram showing the break-points for latitude that we have chosen is shown below ('Graph2' in the Explorer tree).
[![07g._Lat_groupings_hist_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/07g-lat-groupings-hist-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/07g-lat-groupings-hist-i.png)
8. You will want to tweak this new factor for Latitude (just as we did for the Depth factor before). From the NE_Pacific_env data sheet, click **Edit** > **Factors...**, then proceed to do the following:
- ***Duplicate the factor*** of 'Latitude bin ouput' to get 'Latitude bin ouput1' ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07a-duplicate-button.png)).
- ***Rename the factor*** 'Latitude bin ouput1' to call it 'Latitude' ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07b-rename-button.png)).
- ***Rename the levels of the factor*** 'Latitude' to whole numbers ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/07d-rename-levels-button.png)). An image of this last operation is shown below.
[![07f._Rename_Latitude_factor_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/07f-rename-latitude-factor-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/07f-rename-latitude-factor-i.png)
#### Open the groundfish data and import the new factors
Now that we have the factors we want, it would be great to use this to our advantage in analyses of the groundfish data. First we will get the fish data into the workspace, then we will import the new factors (currently associated wtih the environmental data sheet) over to the groundfish data sheet.
9. In the same PRIMER workspace, click **File** > **Open...** and open the file called 'NE_Pacific_groundfish.pri'. It will look like this:
[![08._Groundfish_data_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/08-groundfish-data-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/08-groundfish-data-i.png)
10. From the 'NE_Pacific_groundfish' data sheet in the Explorer tree, click **Edit** > **Factors...** and click the 'Import' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/09a-import-button.png)). Choose to import from the 'NE_Pacific_env' worksheet, then click the 'Select' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/10f-select-button.png)). In the selection dialog, pick only the two factors named 'Depth' and 'Latitude' to include in the import, then click '**OK**', like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/10e-import-operations-all-new.png)
You can confirm that the groundfish sheet now has the 'Depth' and 'Latitude' factors, imported from the environmental data sheet.
#### Create a combined factor of depth-by-latitude
11. For our analysis and plots, we will want now to create a factor that corresponds to the ***combination*** of all depth-by-latitude bins. From the 'NE_Pacific_groundfish', click **Edit** > **Factors...**, then click the 'Combine' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/11b-combine-button.png)). In the 'Combine Factors' dialog, click on the 'Factors...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/12c-factors-button.png)), then in the 'Ordered Selection' dialog, choose to include just Latitude and Depth, as shown below, then click '**OK**' (3 separate times for the three windows).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/12d-combine-factors-all.png)
We now have a factor that identifies groups of samples with similar latitude and depth. This new combined factor of 'Latitude-Depth' effectively corresponds to spatial 'cells' of practical interest in our study design, and it will serve us very well for subsequent analyses. (It is the final column on the right in the image below):
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/12e-finished-all-factor-operations.png)
#### Analyse changes in groundfish assemblages *vs* depth and latitude
We now wish to analyse potential changes in groundfish assemblages with shifts in depth and latitude. Our plan will be to apply a suitable pre-treatment to the data, calculate averages wtihin each latitude-by-depth cell, proceed with calculating a square-root transformation of the data and Bray-Curtis resemblances among these cells, followed by an ordination (nMDS) and tests of hypotheses (ANOSIM) on that resemblance matrix. Note that the averaging step is really important here. We would have no hope of seeing any sensible patterns if we were to 'throw' the full dataset of over 6000 replicates into a single nMDS plot!
### Pre-treatment
To analyse the groundfish data, it is appropriate to consider that many fish species occur in clusters or aggregations of individuals. Hence, a pre-treatment option such as ***dispersion weighting*** ({{2954#bkmrk-clarkeetal2006a}}) would likely be a really appropriate option here. We shall consider the replicate trawls within each latitude-by-depth cell as fairly natural groupings to use in order to apply this pre-treatment option. We noted that there were only 8 replicate trawls from depths less than 50 m, so we will omit those replicates in what follows.
12. From the 'NE_Pacific_groundfish' data sheet, click **Select** > **Samples...** > ($\bullet$Factor levels > Depth), click the 'Levels' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/13d-levels-button.png)), then choose to retain all depth groups except '24', and click '**OK**', as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/13c-select-samples-gt-50m-all.png)
This will turn the worksheet cells blue, and this indicates that a subset of the data has been selected.
[![13e._Subset_selected_groundfish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/13e-subset-selected-groundfish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/13e-subset-selected-groundfish-i.png)
13. From this subset-selected groundfish data sheet, click **Pre-treatment** > **Dispersion Weighting...** and choose to do this on the basis of the factor 'Latitude-Depth', then click '**OK**', like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/14-disp-weighting-dialog.png)
This operation will take some considerable time, simply because of the sheer size of the dataset. The randomization test (for which the individuals of each species are randomly re-assigned to replicates wtihin each of the latitude-depth cells), which is done indpendently for every species, is computationally demanding. However, the results file from the dispersion-weighting pre-treatment (called 'Dispersion weighting1') demonstrates very clearly that ***many*** of these fish species show significant clustering (i.e., wherever the value in the 'Divisor' column is greater than 1), hence should sensibly be pre-treated in this way. The dispersion-weighted data is called 'Data3' in the Explorer tree and will look like this:
[![13_add-on_Dispersion-weighted_data_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/13-add-on-dispersion-weighted-data-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/13-add-on-dispersion-weighted-data-i.png)
### Averaging
14. From the dispersion-weighted data ('Data3'), click **Tools** > **Average...** and choose to average the samples by the factor of 'Latitude-Depth', like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/15-average-dialog.png)
In the resulting data sheet (called 'Data4'), we now have average values for each species in each latitude-by-depth cell (rows), as shown below:
[![15b._Averaged_data_sheet_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15b-averaged-data-sheet-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15b-averaged-data-sheet-i.png)
Note that the names of the samples in the averaged data now correspond to the latitude-by-depth bin combinations. If you click **Edit** > **Properties...**, you will see that this sheet now has 71 rows. We have consolidated these data in a very useful way across our study design, while maintaining the integrity of the underlying information.
### Transformation & resemblance
If we look at a shade plot of the averaged data (you can do this by clicking on **Plots** > **Shade Plot...** from 'Data4'), we can see that, even after dispersion-weighting and averaging, these data still look pretty sparse. We will therefore do a (mild) overall transformation to square roots, then calculate the Bray-Curtis resemblance measure.
15. From 'Data4', click **Pre-treatment** > **Transform(overall)...** > Transformation: Square root, **OK**.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/15c-sqrt-transformation.png)
The resulting square-root transformed data matrx will be called 'Data5'.
16. From 'Data5', click **Analyse** > **Resemblance...** > (Measure $\bullet$ Bray-Curtis similarity) & (Analyse between $\bullet$Samples), '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/15d-bc-resem.png)
The resulting resemblance matrix will be called 'Resem1', and will look like this:
[![15e._Resem_groundfish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15e-resem-groundfish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15e-resem-groundfish-i.png)
### Ordination *via* nMDS
17. From the 'Resem1' matrix, click **Analyse** > **MDS** > **Non-metric MDS (nMDS)...**, retain all of the defaults and click '**OK**'. The resulting best 2d solution for the nMDS ordination plot is called 'Graph3' in the Explorer tree, and it has a nice low stress of 0.06.
We can consider looking at two different 'views' of this ordination:
- **(a)** With symbols corresponding to 'Depth' and Labels correspond to 'Latitude' (optionally with trajectories for 'Latitude', split by 'Depth' groups); or
- **(b)** With symbols corresponding to 'Latitude' and Labels correspond to 'Depth' (optionally with trajectories for 'Depth', split by 'Latitude' groups).
**To obtain (a):**
From 'Graph3', click **Graph** > **Sample Labels & Symbols...** (Labels > $\checkmark$Plot > $\checkmark$By factor Latitude) & (Symbols > $\checkmark$Plot > $\checkmark$By factor Depth). Get the trajectories by clicking **Graph** > **Special**, click the 'Overlays' tab then choose: Trajectory > $\checkmark$Overlay trajectory > Trajectory numeric factor: Latitude > $\checkmark$Split trajectory Depth.
The result looks like this:
[![16a._nMDS_Depth_as_symbols_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/16a-nmds-depth-as-symbols-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/16a-nmds-depth-as-symbols-i.png)
**To obtain (b):**
Simply swap the role of 'Latitude' and 'Depth' factors in the above instructions for (a). The result looks like this:
[![16b._nMDS_Lat_as_symbols_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/16b-nmds-lat-as-symbols-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/16b-nmds-lat-as-symbols-i.png)
These ordinations show a very highly spatially structured ecological system, with marked gradual changes in fish assemblages with both latitude and depth. It is clear that latitudinal turnover in fish asesmblages is more marked at shallower depths than at deeper depths. In addition, turnover in fish assemblages with depth becomes less marked after about 600 m, particularly at higher latitudes.
### Testing ordered factors *via* ANOSIM
We can treat each of these factors as ***ordered factors*** in a ***two-way ANOSIM*** (see {{@954#bkmrk-somerfieldetal2021a}} and {{@954#bkmrk-somerfieldetal2021b}}) to test and quantify these effects in a non-parametric (rank-resemblance) framework.
18. From the 'Resem1' matrix, click **Analyse** > **ANOSIM...**, and specify the two factors of 'Depth' and 'Latitude' as ordered factors in a two-way crossed design, as shown in the dialog below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/17-anosim-dialog-groundfish.png)
The results (not suprisingly) are uber clear (see the output file 'ANOSIM1' in the Explorer tree). There is a highly significant ordered effect of depth generating sequential turnover in groundfish assemblages from 50 m to 1200 m on the NE Pacific coast (ANOSIM, $R^O$ = 0.934, $P$ = 0.0001). There are also significant sequential changes (i.e., turnover) in groundfish assemblages along the latitudinal gradient, from 32°N to 48°N (ANOSIM, $R^O$ = 0.768, $P$ = 0.0001).
[![17b._ANOSIM_output_groundfish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/17b-anosim-output-groundfish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/17b-anosim-output-groundfish-i.png)
[![17c._ANOSIM_Rperm_Depth_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/17c-anosim-rperm-depth-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/17c-anosim-rperm-depth-i.png)
[![17d._ANOSIM_Rperm_Latitude_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/17d-anosim-rperm-latitude-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/17d-anosim-rperm-latitude-i.png)
---
¶*NOAA Fisheries, NWFSC/FRAM, 2725 Montlake Blvd. East, Seattle, WA 98112, USA*
# 15. Other new tools & utilities
# 15.1 New default colour palette
#### Accessibility
It is important to make graphics accessible to those with color vision deficiencies. We have therefore re-vamped the colour palette for P8 to achieve distinctive colours for plots (by default) that carefully accommodate the most common forms of colour blindness (Fig. 15.1).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/02b-colour-palette-comparison.png)
*Fig. 15.1 New default colour palette in PRIMER 8 (at right) compared to the old default colour palette in PRIMER 7 (at left).*
In developing these new defaults, we used the following resources:
- ["Coloring for Colorblindness" - David Nichols](https://davidmathlogic.com/colorblind/#%23D81B60-%231E88E5-%23FFC107-%23004D40)
- ["Points of view: Color blindness" - Bang Wong](https://www.nature.com/articles/nmeth.1618)
- [Paul Tol's Notes: Colour schemes and templates - Paul Tol](https://sronpersonalpages.nl/~pault/)
We especially appreciated David Nichols' [online tool](https://davidmathlogic.com/colorblind/#%23D81B60-%231E88E5-%23FFC107-%23004D40), which allowed us to play with colour schemes and get a feel for how different schemes would look for people who experience colour differently, due to various forms of alterations to their color vision, including protanopia, deuteranopia, tritanopia or deuteranomaly. We ended up with the following pallette:
- [PRIMER 8 default colours](https://davidmathlogic.com/colorblind/#%230072B2-%23E69F00-%236ED065-%2356B4E9-%23F0E442-%23D55E00-%23CC79A7-%23009E73-%239c0ed3-%23934f38-%23c8e761) - (see the left-hand column labeled 'True' when you follow this link).
To get there, we borrowed the base colour scheme from {{@954#bkmrk-wong2011}}. We then simply re-arranged the order of the colours and extended the scheme to 12 colours which we felt, as much as possible, would still appear well-differentiated from one another, right across the range of different forms of colour vision deficiency.
#### Compatibility with earlier versions of PRIMER
If you have an existing *.pri or *.pwk file that has been opened and saved in PRIMER 7, then it will carry around the default colour palette from PRIMER 7 with it. Thus, when you open up a PRIMER 7 file in PRIMER 8, those old default colours from PRIMER 7 (or whatever other colours you tweaked and saved in your file's graphics) will be used. To get rid of all of the older colours, and hence see the new default colour palette in P8 when you create a new graphic, please save the file as an Excel file and then import the data, fresh, into PRIMER 8.
# 15.2 New selection options
In PRIMER 8, the options available for ***selecting samples*** or ***selecting variables*** have been expanded considerably from what they were in PRIMER 7.
#### Selecting Samples
When you click **Select** > **Samples...** in P8, you will see the new dialog shown at right in Fig. 15.2.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/03c-select-samples-compare.png)
*Fig. 15.2 Comparison of the **Select** > **Samples...** dialog window in PRIMER 7 (at left) vs PRIMER 8 (at right).*
Note that in PRIMER 8, you can select samples:
- using **sample names** (this option includes a filter to help find names in long lists);
- using **sample numbers**;
- belonging to certain **level(s) of a factor**; or
- that contain some minimum number ($x$) of **non-zero values**.
You can also ***exclude*** samples that:
- have ($x$) or more zero values;
- have ($x$) or more missing values; or
- have ($x$)% or more missing values.
Note also that there is the option to:
'$\checkmark$Output selection to new worksheet'.
This latter option means you don't have to take any extra steps (post-selection) simply to duplicate the selected data into a new worksheet (if desired) for subsequent analyses. As an added bonus, a small output file is produced when you tick this option which identifies precisely what choices you made to select the samples. It is very useful to have this information going forward, to keep track of any subset selections made along your analysis pathway.
#### Selecting Variables
When you click **Select** > **Variables...** in P8, you will see the new dialog shown at right in Fig. 15.3.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/04c-select-vars-compare.png)
*Fig. 15.3 Comparison of the **Select** > **Variables...** dialog window in PRIMER 7 (at left) vs PRIMER 8 (at right).*
Note that in PRIMER 8, you can select variables:
- using **variable names** (this option includes a filter to help find names in long lists);
- using **variable numbers**; or
- belonging to certain **group(s) of an indicator**.
It is also possible to select variables that correspond to the top ($x$) variables in a list based on:
- **frequency of occurrence**;
- **total abundance** (sum); or
- **percent contribution to any one sample** (this option was called 'most important' in PRIMER 7)
Finally, you may choose to select variables that:
- occur in at least ($x$) sample(s);
- occur in at least ($x$)% of sample(s); or
- contribute at least ($x$)% to the total abundance (sum).
Note that the above dialog gives you a lot of freedom with respect to identification of 'important' variables in ecological contexts, not just on the basis of total abundance (in any one sample or overall), but alternatively by reference to their frequency of occurrence. Using the frequency of occurrence may often be more suitable, as a criterion, than pulling out species based on their total abundance values, particularly if there are highly sporadic species that have massive abundance values (e.g., weeds or opportunists), but that may not actually be that important ecologically (e.g., they might have shown up in just a single sample).
As an added bonus (and precisely as we saw in the dialog for the selection of samples, shown above), you can choose to:
'$\checkmark$Output selection to new worksheet'.
This streamlines and clarifies analytical pathways.
#### Filters *(New!)*
For either the selection of samples or the selection of variables **by name**, PRIMER 8 offers a new tool in the form of a ***filter*** that is very helpful whenever you are dealing with long lists of names in large data sheets. For example, suppose I wish to select all of the variables from a long list of species that belong to the same genus (e.g., *Ampelisca*), I simply choose to select variables by names and click the 'Select Variable Names...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/05d-select-var-names-button.png)), then start typing the genus name into the 'Filter:' box, and those species will appear in the 'Available:' list for me to easily select them (see Fig. 15.4 below).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/05c-select-vars-with-filter.png)
*Fig. 15.4 Dialog window for selecting names of variables, using the new 'Filter:' tool in PRIMER 8. In the above example (data in the file named 'Norway_macrofauna.pri', found inside the 'Examples_P8 > Norway_macrofauna' folder), there are 809 taxa in the data sheet, but by typing the letters 'Ampe' in the 'Filter:' box, the available list is reduced to just the 12 variables that have this particular combination of letters in their name. All of these are of the genus 'Ampelisca'.*
Using this new filtering tool in PRIMER 8 means you can quickly and easily find and select specific variables (or samples) that you know are in your data somewhere (e.g., all variables starting wtih the letter "R"), even if you cannot easily remember the entire name in detail.
# 15.3 Re-name levels of a factor (or indicator)
There are many situations where it would be very handy to be able to change the names of levels of a factor (or to change the names of groups for an indicator). We often need to tweak the names of levels of factors. For example,
- The names are too long and you want to shorten/abbreviate them so that they are easy to see as labels on an ordination plot. For example, you have levels of 'Berghan Point' and 'Home Point', and you want to change them to 'B' and 'H'.
- The factor has levels that are quantitative (e.g., such as depth zones), but the names are not strictly numeric (e.g., maybe the levels are named '50-100 m', '100-200 m', '200-300 m', etc.). You might, however, really want numeric factor level names so you can analyse the factor as **ordered** in an ANOSIM, put depth trajectories on ordination plots, etc.
- You have discovered that not all factor levels matter, and so you want to combine 2 or more factor levels into an amalgamated (single) factor level. For example, maybe you originally have levels of 'Low', 'Medimum-Low', 'Medium', and 'High', but later decide that the first two levels are not actually significantly different, so you want to combine those two and change the names so that you just have 'L', 'M' and 'H'.
- You have taken data through time (say, every 2 months) at your sites, and your labels for this temporal factor look like this: 'Jan-2001', 'Mar-2001', 'May-2001', etc., but you want to treat these simply as sequential time-points in your plots, e.g., 'T1', T2', 'T3', etc.
It would be very laborious to have to go through the process of re-naming the factor-level names for every single sample in the full worksheet, and this process would also be highly prone to error.
In PRIMER 8, all you need to do is click **Edit** > **Factors**, click on the factor whose level names you want to change, and then click the 'Rename Levels...' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/06d-rename-levels-button.png)). The resulting dialog lets you nominate the new names for the factor levels (just type them in), and then you can get on with your work.
For example, for the size-class data for mussels from several sites in the Gulf of Alaska (the data file is 'Gulf_of_Alaska_mussels.pri', located in the folder 'Examples_P8' > 'Gulf_of_Alaska_mussels'), we may wish to abbreviate the names of the ecoregions so that they are shorter, making them easier to see in ordination plots (Fig. 15.5).
[![06c._Factors_rename_levels_mussels_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/06c-factors-rename-levels-mussels-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/06c-factors-rename-levels-mussels-i.png)
*Fig. 15.5 Dialog windows showing the action of changing the names of levels for the factor of 'Ecoregions' for the Gulf of Alaska mussel size-class dataset: from 'Lynn Canal' to 'LC' and from 'Kachemak Bay' to 'KB'.*
The resulting factor information (after changing the level names for 'Ecoregion', as shown above) looks like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/06e-after-changing-level-names.png)
As an aside, if you wanted to retain the original level names in your data sheet as well, then you need first to ***duplicate*** the original factor (using the 'Duplicate' button), ***rename*** the factor itself (using the 'Rename...' button; for example, the abbreviated names could be held in a factor called 'Ecoreg'), ***then*** create new level names for that duplicated factor (using the 'Rename Levels...' button) to whatever new names you wish.
You can compare the metric MDS plot of the sites (based on Manhattan distances among cumulative percentages of standardised size-classes of mussels) when Ecoregion factor-level labels have long names *vs* the same plot with short names (see below).
[![06f._mMDS_mussels_full_names_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/06f-mmds-mussels-full-names-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/06f-mmds-mussels-full-names-i.png)
[![06g._MDS_mussels_short_names_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/06g-mds-mussels-short-names-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/06g-mds-mussels-short-names-i.png)
Finally, note that this same re-naming tool is also available under **Edit** > **Indicators...**.
# 15.4 Add customised values/labels to graphical axes
In PRIMER 8, there is a new tool that allows us to ***add customised values and labels*** to coordinate axes in graphics.
Consider the following scatter plot of species richness ($S$) *vs.* structural complexity of the substratum (measured using a chain-and-tape method), obtained from a visual survey of intertidal macrofauna and algae at a site inside the marine reserve at Long Bay, on the north shore of Auckland, New Zealand.¶ Biotic data ('Long_Bay_intertidal_biota') and environmental data ('Long_Bay_intertidal_env.pri') from this study are found in the 'Examples_P8' > 'Long_Bay_intertidal' folder.
[![07a._Long_Bay_Scatter_S_vs._Complex_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/07a-long-bay-scatter-s-vs-complex-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/07a-long-bay-scatter-s-vs-complex-i.png)
In the above scatter plot, the points that correspond to the minimum ('Min') and maximum ('Max') values obtained for the variable of 'Complexity' (on the x-axis), have been labeled individual (using a factor that is empty for all other samples). We might decide it would be useful to annotate the graphic to provide the minimum and maximum values on the x-axis for those two values directly as well. We can do this by clicking on the x-axis (or by clicking **Graph** > **General...** and clicking on the 'X-axis' tab).
In this dialog, we can see an 'Additional Labels' button ([](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/7c-additional-labels-button.png)), which is **new in P8**. We can also see the new option to change the **label orientation** so that it is either '$\bullet$Perpendicular' or '$\bullet$Parallel'. For the present example, let's choose '$\bullet$Parallel' and then go ahead and click the 'Additional Labels' button where we can add both the specific **values** and the desired **labels** for the maximum and minimum complexity, also including a tick-mark for each of them, like so:
[![7e._Additional_Labels_ALL_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/7e-additional-labels-all-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/7e-additional-labels-all-i.png)
The resulting graphic looks like this:
[![07f._Long_Bay_Scatter_with_Additional_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/07f-long-bay-scatter-with-additional-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/07f-long-bay-scatter-with-additional-i.png)
There are clearly many uses for this new feature. It is a great way to provide clarity regarding any specific values (or sets of values) that you may wish to highlight along particular axes in customised graphics.
### Change the axis label orientation
Note that another new feature in PRIMER 8 is the possiblity to change the ***orientation*** of the axis labels. See [section 3.2](https://learninghub.primer-e.com/link/1023#bkmrk-change-the-axis-labe) for an example.
---
¶*These data were collected in 2015 as part of a course in quantitative marine ecology at Massey University, Albany, New Zealand. Taxa that were recorded in the survey as 'dead' (e.g., oyster shells, barnacle tests, etc.) were not included in the tally of richness analysed here.*
# 15.5 Split data sheet by factor/indicator
In PRIMER 8 there is a new tool, accessed by clicking **Tools** > **Split Data...**, which allows you to split a data sheet into several separate data sheets, corresponding to:
- separate groups of ***samples***, based on a ***factor***; or
- separate groups of ***variables***, based on an ***indicator***.
This is a simple tool, but a helpful one. Having this tool saves one from having to select each group individually, then duplicate them into a new sheet, one group at a time.
For example, consider the data set of fish surveys from New Zealand, discussed in [section 10.4](https://learninghub.primer-e.com/link/1049) above. (Data are located in the file 'NE_NZ_fish_counts.pri', found in the 'Example_P8' > 'NE_NZ_fish' folder.) Visual underwater surveys of fish assemblages were done at four different locations along the north-eastern coast of New Zealand: Berghan Point, Home Point, Leigh and Hahei. Suppose now we would like to examine trends in assemblages through time, but do this separately for each of these four different locations.
From the fish dataset, click **Tools** > **Split Data...**, like so:
[![08a._Tools_Split_Data_fish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/08a-tools-split-data-fish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/08a-tools-split-data-fish-i.png)
We can then choose to split the data based on the factor of 'Location', as shown in the 'Split' dialog window below, then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/8b-split-data-dialog-window-fish.png)
In our Explorer tree, we then see four new data sheets, one for each of the four locations:
[![08c._Post_data_split_fish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/08c-post-data-split-fish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/08c-post-data-split-fish-i.png)
Note that the sheets produced by this operation each have a new ***title*** (see also **Edit** > **Properties**) that reflects:
- the factor on which the split was made; and
- the relevant corresponding group for each resulting data sheet.
For example, the title for 'Data4' is given as 'New Zealand Fish Visual Surveys - Location: Home.Point', *viz.:*
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/08d-new-title-home-point.png)
Of course, you can also change the name of the data sheets themselves if you wish. For example, you can change the name of the sheet called 'Data4' to 'Home.Point' by clicking **File** > **Rename Data**.
# 15.6 Line plots for samples
There is a new facility in PRIMER 8 to create **Line plots** in two different ways:
- with one line for every ***variable*** (across all samples); or
- with one line for every ***sample*** (across all variables).
This is a considerable improvement on the line plot dialog in PRIMER 7, which only could be implemented to draw lines for variables (Fig. 15.6).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/09e-line-plot-compare-p7-p8.png)
*Fig. 15.6 Comparison of the **Plots** > **Line Plot...** dialog window in PRIMER 7 (at left) vs PRIMER 8 (at right).*
We saw this tool earlier in [section 13.3](https://learninghub.primer-e.com/link/1063), in the analysis of mussel size-classes at different sites in the Gulf of Alaska (the data file is called 'Gulf_of_Alaska_mussels.pri', located in the folder 'Examples_P8' > 'Gulf_of_Alaska_mussels').
From the mussel data sheet, after standardising the original raw count data to [cumulative percentages](https://learninghub.primer-e.com/link/1063#bkmrk-input-data-and-stand) (the variables are the sizes classes here), resulting in a sheet called 'Data1', we can choose to create a line plot with a line for each sample (the sites here) by clicking **Plots** > **Line Plot...**, like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/09c-line-plot-menu-item.png)
In the 'Line Plot'dialog window, we simply choose to view lines for the ($\bullet$ Samples), and then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/09b-line-plot-p8.png)
*Note: change the above dialog if the wording changes*
Note in the above dialog that we also have the option to output ***multiple line plots*** corresponding to groupsof samples identified *via* a factor (or groups of variables identified by an indicator, if we are drawing lines for variables).
The resulting line plot for this example is shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/09d-line-plot-mussel-data-take2.png)
# 15.7 Output group-level stats from dispersion (or variability) weighting
In PRIMER 8, you can now output more detailed statistical information from either dispersion-weighting or variability weighting to a worksheet. You can see this additional option in the comparison of the 'Dispersion Weighting' dialog in PRIMER 8, compared to that in PRIMER 7 (Fig. 15.7).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/10c-dw-comparison-p7-p8.png)
*Fig. 15.7 Comparison of the **Pre-treatment** > **Dispersion Weighting...** dialog window in PRIMER 7 (at left) vs PRIMER 8 (at right).*
In essence, whenever we calculate the average ***index of dispersion*** ($\bar{D}$) for the dispersion weighting pre-treatment, it is drawn from individual indices of dispersion (variance-to-mean ratios, $D_i$) for each of several ($i = 1, \ldots, g$) groups (identified by a factor). It would be helpful, in some situations, to be able to see those original individual indices of dispersion across all of the groups for each variable. The extra tickbox ('$\checkmark$Output group-level worksheet') available in the PRIMER 8 dispersion weighting dialog allows you to produce and examine this underlying statistical information in a worksheet directly.
For example, we can re-visit the data we saw in [section 10.4](https://learninghub.primer-e.com/link/1049) above, located in the file 'NE_NZ_fish_counts.pri', found in the 'Example_P8' > 'NE_NZ_fish' folder. We can begin by following steps 1, 2 and 3 from [section 10.4](https://learninghub.primer-e.com/link/1049) on these data; i.e., get the data into PRIMER 8, then:
- Obtain a subset of the data based on the factor of 'Year', from 2010-2015 inclusive (**Select** > **Samples...**), with the resulting subset being called 'Data1'.
- Create a combined factor consisting of all combinations of the factors: 'Loc', 'Hab' and 'Year' (**Edit** > **Factors...** > **Combine...**)
Next, we can apply the dispersion weighting pre-treatment on the basis of this combined factor by clicking **Pre-treatment** > **Dispersion Weighting...** and in the resulting dialog, choose Factor: Loc-Hab-Year. This time, however, we shall also choose to get ($\checkmark$Stats to worksheet) & ($\checkmark$Output group-level worksheet), then click '**OK**'.
This will produce the following three worksheets (assuming 'Data1' contains the subset-selected data):
- 'Data2' contains the dispersion-weighted data.
- 'Data3' contains the statistics associated with the ***average*** index of dispersion ($\bar{D}$) for each variable.
- 'Data4' contains the ***individual*** indices of dispersion ($D_i$) for each group (columns) and for each variable (rows).
These three data sheets are shown below for the present example.
[![10d._DW_data_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/10d-dw-data-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/10d-dw-data-i.png)
[![10e._D-bar_values_fish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/10e-d-bar-values-fish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/10e-d-bar-values-fish-i.png)
[![10f._D_i_values_fish_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/10f-d-i-values-fish-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/10f-d-i-values-fish-i.png)
Note in the above sheet ('Data4', containing individual indices, $D_i$) that there are many missing values. This is simply a consequence of there being zero fish of that species in that particular group of samples; hence, the mean and the variance will both be equal to zero, so no value of $D_i$ can be calculated for that particular cell (or group). Average $\bar{D}$ values for each species are naturally calculated only from those cells (groups of samples) where at least some individuals of that species were recorded (i.e., the non-missing entries).
Finally, note that this '$\checkmark$Output group-level worksheet' option is also available for **Pre-treatment** > **Variability Weighting...**, *viz.:*
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/10g-variability-weighting-dialog-too.png)
# 15.8 Output diagnostic plots from CAP
In PRIMER 8, it is now possible to output diagnostic plots from the CAP routine (i.e., canonical analysis of principal coordinates, an analysis obtained from a resemblance matrix by clicking **PERMANOVA+** > **CAP**). You simply tick the new option to '$\checkmark$Do diagnostic plots' in the CAP dialog window.
For example, consider a study of the assemblages of small benthic fishes living in rocky subtidal reef habitats in northeastern New Zealand, examined by {{@954#bkmrk-smithanderson2016}}. Data from this study are located in the file called 'NZ_benthic_fish.pri', in the 'Examples_P8' > 'NZ_benthic_fish' folder. Surveys were done of the benthic fish fauna, along with fine-scale habitat features in kelp forests and rocky reefs along the north-eastern coast of New Zealand. There were sites at a range of locations in and around several marine reserves (including Leigh, Tawharanui, Hahei and the Poor Knights Islands), and data were obtained over a period of 3 years (2011-2013). At each site, divers surveyed $n$ = 8 transects, measuring 1 m $\times$ 5 m. Each transect was made up of five contiguous 1 m $\times$ 1 m quadrats. For each quadrat, divers first visually searched and recorded counts for all benthic fishes, then recorded the presence/absence of a set of pre-defined habitat features (see the file named 'NZ_benthic_fish_habitat.pri', located in the same folder).
Here, we shall consider only the potential effects of the marine reserve at Leigh on these fish communities in a fully balanced design. Note that the factor of 'Reserve' has two levels: 0 = outside the reserve, and 1 = inside the reserve. Due to the sparsity of the count data at small spatial scales, we shall analyse data summed to the transect level, then consider densities of each species (averages per transect) at the spatial scale of sites for the ensuing CAP analysis. Mean densities (per 5m2) of each fish species per site are located in the file named 'NZ_benthic_fish_Leigh_av_densities.pri'.
1. Open the file ('NZ_benthic_fish_Leigh_av_densities.pri') in PRIMER 8. It will look like this:
[![12a._Leigh_fish_densities_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/12a-leigh-fish-densities-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/12a-leigh-fish-densities-i.png)
2. From 'NZ_benthic_fish_Leigh_av_densities', click **Analyse** > **Resemblance...** and choose to calculate the Bray-Curtis similarity among samples. This will produce a resemblance matrix called 'Resem1'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/12b-resem-tfins.png)
[![12bb._Resem_matrix_tfins_[i]_NEW.png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/12bb-resem-matrix-tfins-i-new.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/12bb-resem-matrix-tfins-i-new.png)
3. From the 'Resem1' matrix, run the CAP routine aiming to distinguish small benthic fish communities inside *vs* outside the marine reserve, by clicking **PERMANOVA+** > **CAP** and choose:
(Analyse against $\bullet$Groups in factor > Factor for groups or new samples: 'Reserve') &
($\checkmark$Scores to worksheet) &
(Diagnostics > $\checkmark$Do diagnostics > $\checkmark$Do diagnostic plots) &
($\checkmark$Do permutation test > Num. permutations: 9999)
then click **OK**, as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-07/12c-cap-dialog-tfins-p8.png)
The option to ***do diagnostic plots*** is new to PRIMER 8.
The analysis will run and produce the following:
- 'CAP1' ([![13._CAP1_notepad_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/13-cap1-notepad-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/13-cap1-notepad-i.png)), containing all of the essential results in a rich text format (*.rtf), and
- 'MultiPlot1' ([![14._MultiPlot1_Exp_tree_icon_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/14-multiplot1-exp-tree-icon-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/14-multiplot1-exp-tree-icon-i.png)), containing the CAP plot itself (as 'Graph1') along with plots of the ***diagnostics for choosing $m$***, where $m$ = the number of principal coordinate axes (PCOs) used to produce the canonical (CAP) axis (as 'Graph2' through 'Graph5').
From 'CAP1' and the accompanying CAP plot ('Graph1'), shown below, we can see that the canonical correlation associated with the 'Reserve' effect is reasonably strong ($\delta_1$ = 0.644), and that this is statistically significant ($\delta_1^2$ = 0.415, $P$ = 0.001). This CAP model was achieved with $m$ = 4 PCO axes, which achieved a leave-one-out allocation success of 77.78% under cross-validation.
[![12d._CAP_rtf_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/12d-cap-rtf-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/12d-cap-rtf-tfins-i.png)
[![15a._Graph1_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15a-graph1-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15a-graph1-tfins-i.png)
The four ***diagnostic plots*** are shown below.
[![15b._Graph2_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15b-graph2-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15b-graph2-tfins-i.png)
[![15c._Graph3_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15c-graph3-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15c-graph3-tfins-i.png)
[![15d._Graph4_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15d-graph4-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15d-graph4-tfins-i.png)
[![15e._Graph5_tfins_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/15e-graph5-tfins-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/15e-graph5-tfins-i.png)
It is clear that the choice of $m$ = 4 is a good one here, as this choice achieved the greatest leave-one-out allocation success ('Graph2') and the lowest leave-one-out residual sum-of-squares ('Graph3'). Although this can also (technically) be seen in the '*DIAGNOSTICS*' section of the 'CAP1' output above, it is certainly much clearer to see this information in diagnostic plots, whereby the wisdom of the choice made, overall, by reference to other values of $m$, can be readily assessed.
# 15.9 New diagnostics for PCA/PCO plots
#### Background
Consider a cloud of $N$ points (sampling units) in a $p$-dimensional multivariate space. [Principal components analysis (PCA)](https://learninghub.primer-e.com/link/102) is an ordination method that will find an axis through that cloud of points so as to maximise the total variance of the points, when they are projected at right angles onto that axis. Having found such an axis, a second axis is then built in a similar way, subject to it being perpendicular to (independent of) the first axis, and so on. PCA does this in Euclidean space. [Principal coordinate analysis (abbreviated as PCO or PCoA)](https://learninghub.primer-e.com/books/permanova-for-primer-guide-to-software-and-statistical-methods/chapter/chapter-3-principal-coordinates-analysis-pco) also does this, but in the space of a chosen resemblance measure ({{@954#bkmrk-gower1966}}).
Hence, we may conceptually consider either PCA or PCO as a type of ***projection*** of the points into a smaller number of dimensions that can capture important major stuctures occurring across the data cloud as a whole. Specifically, if there is some redundancy of information (inter-correlation) among the $p$ original variables, then we may anticipate that a subset of (say) two or three principal axes will capture a substantial portion of the total variation in the system. A configuration plot of the positions of the points along the first two (or three) principal axes can therefore be examined in an attempt to visualise patterns among the sampling units along these major axes of variation.
To assess the utility of a 2-d (or 3-d) PCA/PCO plot, one might simply examine the percentage of the total variation that is extracted by the first 2 (or 3) axes. If this is substantial (e.g., 70% or 80%), then we may consider the plot to be useful for interpretation. Another useful tool is a [***scree plot***](https://learninghub.primer-e.com/link/753), which shows the percentage of variation explained by sequentially ordered principal axes.
This does not, however, provide us with a meaningful measure of how well the plot (of reduced dimension) represents the (high-dimensional) inter-point distances (or dissimilarities) in the original multivariate space.
A good way to do that would be to calculate the [***stress***](https://learninghub.primer-e.com/link/107#bkmrk-measure-goodness-of-) associated with the 2-d (or 3-d) configuration. The notion of stress is already quite familiar as a measure of the adequacy of an ordination plot produced using multi-dimensional scaling (MDS, see [section 5.2](https://learninghub.primer-e.com/link/107) and [section 5.8](https://learninghub.primer-e.com/link/114) in *[Change in Marine Communities, 3rd edition](https://learninghub.primer-e.com/books/change-in-marine-communities)*). Indeed, it is the specific goal of the MDS algorithm to create a configuration that minimises stress. Of course, neither PCA nor PCO are specifically designed to minimise stress, so a PCA or PCO configuration will always (necessarily) have higher stress than an MDS configuration for a given data set. Nevertheless, having a measure of stress for a PCA/PCO configuration can help us assess its utility for displaying inter-point relationships faithfully. Here, we would consider that the usual rule-of-thumb of stress < 0.20 indicates a usefully interpretable configuration. Also, a plot of the configuration distances *vs* the original inter-point distances - [a Shepard diagram](https://learninghub.primer-e.com/link/107#bkmrk-specify-the-number-o) - is clearly also a useful diagnostic tool here.
#### Scree Plot, Shepard diagram and Stress
In the PCO routine (accessed from a resemblance matrix by clicking on **PERMANOVA+** > **PCO...**), PRIMER 8 now offers you new options to output these diagnostics tools:
- a Scree plot;
- Shepard diagrams in 2D and/or 3D; and
- Stress calculated using either metric MDS or threshold metric MDS
Fig. 15.8 compares the dialog window for the PCO routine in PRIMER 7 versus PRIMER 8.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/16-new-pco-dialog.png)
*Fig. 15.8 Comparison of the **PERMANOVA+** > **PCO...** dialog window in PRIMER 7 (at left) vs PRIMER 8 (at right).*
These same new diagnostic output options are also available for the PCA routine in PRIMER 8 (obtained from a data matrix by clicking on **Analyse** > **PCA...**).
It is important to recognise that these additional diagnostics do not in any way change anything about the fundamental PCA or PCO analysis itself that is being done. The calculation of principal axes for PCA or PCO ordination still remains exactly the same as ever, and may best be thought of as a *projection*. What is new and different here is simply the opportunity also to view the stress and the Shepard diagrams which attend the resulting PCA or PCO configuration. These diagnostics also help put these methods into perspective vis-à-vis any MDS ordination solutions one might obtain for a given dataset.
#### Example: Messolongi lagoon diatoms
Let's consider a study of diatom assemblages (densities of 193 species) at 17 sites in the lagoons of Messolongi, Aitoliko and Kleissova in Eastern Central Greece ({{@954#bkmrk-danielidis1991}}). Data from this study are contained in the file called 'Messolongi_diatom_density.pri', located in the 'Examples_P8' > 'Messolongi_diatoms' folder.
In what follows, we shall create an ordination of these data on the basis of a Bray-Curtis resemblance matrix calculated from square-root transformed densities, using each of the following methods:
- Non-metric MDS (nMDS)
- Threshold-metric MDS (tmMDS)
- Metric MDS (mMDS); and
- Principal co-ordinates analysis (PCO).
We will compare and contrast not only the ordination diagrams produced, but also the associated diagnostic Shepard diagrams and the values of stress. For PCO, we can produce Shepard diagrams and stress values in two different ways: (i) for comparison with tmMDS; or (ii) for comparison with mMDS.¶
1. Open up the data matrix (Messolongi_diatom_density.pri) in PRIMER 8.
[![17._Messolongi_data_only_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/17-messolongi-data-only-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/17-messolongi-data-only-i.png)
2. From the data matrix (Messolongi_diatom_density), click **Pre-treatment** > **Transform(overall)...** and choose Square-root.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/18-sqrt-transform.png)
3. The square-root transformed data are held in the sheet named 'Data1'. From this, click **Analyse** > **Resemblance...** and choose to calculate Bray-Curtis resemblances among samples, like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/19-bc-resem.png)
### Non-metric MDS (nMDS)
4. The Bray-Curtis similarities are held in the sheet called 'Resem1'. From this, click **Analyse** > **MDS** > **Non-metric MDS (nMDS)...**. Change the 'Max. Dimension' from 3 to 2 (for this exercise, we shall focus only on the 2-d solution), then click 'OK', as shown in the dialog below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/20-nmds-dialog.png)
In the results ('MultiPlot1'), we can see the nMDS ordination ('Graph1'), the ***monotonic*** relationship shown in the Shepard diagram ('Graph2'), and that the value of stress is 0.092.
[![21._nMDS_results_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/21-nmds-results-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/21-nmds-results-i.png)
### Threshold metric MDS (tmMDS)
5. From 'Resem1', click **Analyse** > **MDS** > **Metric MDS (mMDS / tmMDS)...**, then choose 'Choice of intercept: $\bullet$Threshold metric MDS (non-zero intercept)', change the 'Max. Dimension' from 3 to 2, then click 'OK' (see below).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/22-tmmds-dialog.png)
In the results ('MultiPlot2'), we can see the tmMDS ordination ('Graph3'), the ***linear*** relationship shown in the Shepard diagram ('Graph4'), with a ***non-zero intercept*** of 60.54% similarity, and a stress value of 0.130. This is a bit larger than the stress obtained for the nMDS (0.092), but is still sufficiently low (< 0.20) to permit useful interpretations of the patterns of inter-point relationships seen in the diagram.
[![23._tmMDS_results_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/23-tmmds-results-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/23-tmmds-results-i.png)
Note that the threshold (intercept) value of 60.54% similarity indicates that two points that are 'on top of one another' in the tmMDS configuration are to be interpreted as being *ca*. 39.46% dissimilar, and not 0% dissimilar.
Also, unlike the nMDS plot (which has no labels on its axes and from which only rank-order relationships among the points can be interpreted), the ***axes are labeled*** on the tmMDS plot, with values given in the units of the original (Bray-Curtis) resemblance measure. Thus, two points that are (say) 20 units apart from one another on this plot can be interpreted as being approximately (20% + 39.46%) = 49.46% dissimilar. Thus, although the tmMDS comes at the price of higher stress compared to the nMDS, it also 'gives something back' in the form of added interpretability (i.e., with labels on the axes permitting the estimation of actual inter-point dissimilarities). So, clearly, tmMDS can be a good option, provided the stress still remains low enough to yield a plot that is adequate for interpretation.
### Metric MDS (mMDS)
6. From 'Resem1', click **Analyse** > **MDS** > **Metric MDS (mMDS / tmMDS)...**, then choose 'Choice of intercept: $\bullet$Metric MDS (zero intercept)', change the 'Max. Dimension' from 3 to 2, then click 'OK', viz:
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/24-mmds-dialog.png)
In the results ('MultiPlot3'), we can see the mMDS ordination ('Graph5'), the ***linear*** relationship with a ***zero intercept*** in the Shepard diagram ('Graph6'), and a stress value of 0.252. This is much larger than that of either the nMDS (0.092) or the tmMDS (0.130), and is so high (> 0.25) as to prevent meaningful interpretation of any patterns seen in the resulting diagram.
[![25._mMDS_results_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/25-mmds-results-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/25-mmds-results-i.png)
### Principal Co-ordinate Analysis (PCO)
7. First, we will run the PCO and get diagnostics for a Shepard diagram and associated stress calculation that uses a zero intercept (as in mMDS). From 'Resem1', click **PERMANOVA+** > **PCO**, then choose 'Shepard diagram > Stress calculated using: $\bullet$mMDS (1-to-1)', then click 'OK' (see the dialog below).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/26-pco-w-mmds-shepard-dialog.png)
The resulting PCO plot ('Graph7') and associated diagnostic Shepard diagram ('Graph8') are shown in 'MultiPlot4' in the Explorer tree.
[![27._PCO_w_mMDS_Shepard_results_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/27-pco-w-mmds-shepard-results-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/27-pco-w-mmds-shepard-results-i.png)
The first two principal axes explain 44.56% of the total variation (less than half, as seen in the output file called 'PCO1' in the Explorer tree). In addition, the stress is extremely high, at 0.30, indicating that many of the inter-point distances seen in the PCO plot do not reflect their true dissimilarity well at all. This is especially true for the small to middle-sized dissimilarities that range from about 40% to 70% (see the large scatter of points in 'Graph8' for values having ~ 60% to 30% similarity along the x-axis).
8. Next, we will run the PCO again, but this time choose 'Shepard diagram > Stress calculated using: $\bullet$tmMDS' (see the dialog below).
[](https://learninghub.primer-e.com/uploads/images/gallery/2025-11/28-pco-w-tmmds-shepard-dialog.png)
The results are given in 'MultiPlot5' in the Explorer tree.
[![29._PCO_w_tmMDS_Shepard_results_[i].png](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/scaled-1680-/29-pco-w-tmmds-shepard-results-i.png)](https://learninghub.primer-e.com/uploads/images/gallery/2025-12/29-pco-w-tmmds-shepard-results-i.png)
It is clear that the PCO plot itself remains exactly the same (compare 'Graph9' with 'Graph7'). However, we are now using a different yardstick to measure the stress associated with this plot. Specifically, we are permitting the intercept to float away from the origin, and asserting (as we did with the tmMDS) that there is a ***threshold dissimilarity*** that a pair of samples must exceed before they will occupy two different positions on the configuration plot. Here, the intercept is a similarity of 63.24%, corresponding to a dissimilarity threshold of 36.76%. In other words, any two samples that appear on top of each other on the plot (at a distance of zero) can be estimated to be 63.24% similar, not 100% similar.
By adding this flexibility, we can see that the stress improves dramatically, dropping to 0.201, which borders on becoming acceptably interpretable, at least for the broader structures. Once again, it is the larger dissimilarities that are captured better by the PCO (see the tighter clustering of points around the line-of-best-fit in the upper right area of the Shepard diagram), compared to the smaller ones.
### Observations & recommendations
This example serves to demonstrate a more general point. Namely, for a given dataset, the ordering of methods with respect to the value of stress associated with the final configuration will almost always be: nMDS < tmMDS < mMDS < PCO. This is so because:
* all MDS solutions will, by their very design, achieve lower stress than a projection-type ordination, such as PCO, as they are explicitly geared to minimise stress; and
* within the MDS family of methods (nMDS, tmMDS and mMDS), nMDS has the fewest constraints, requiring only a monotonic relationship, while mMDS has the greatest constraints, requiring not just a linear relationship but also a zero-intercept in the Shepard diagram. This is more difficult to achieve, resulting in a higher stress.
Our general recommendations for ordination are:
* Use nMDS, as starting point, to achieve the lowest-stress solution.
* If the Shepard diagram for the nMDS shows that the relationship between configuration distances and original distances is approximately linear, then tmMDS can be used to achieve an ordination with the added benefit of having axes which (along with the added threshold dissimilarity) can be used to estimate actual inter-point dissimilarities directly from the plot. Only retain the tmMDS, however, if the stress is sufficiently low to permit interpretability ( < 0.20 as a rule-of-thumb).
* Only go one step further to consider using mMDS if the intercept is sufficiently close to zero to warrant this. (This will be evident in the Shepard diagram).
* PCO will typically ***not*** be the best tool to use, in general, to visualise inter-point relationships in an ordination, but ***will*** extract the axis (or axes) of greatest total variation through the data cloud as a whole. This means that large dissimilarities will typically be represented much better than small dissimilarities in a PCO plot.
### Utility of PCO
Although PCO may not be the tool of choice for ordination, it is important to recognise that it is still a very useful tool for other reasons. Specifically, a variety of PERMANOVA+ routines use PCO ('under the hood') in order to calculate correctly (using the ***full*** set of PCO axes, and carefully accounting for negative eigenvalues):
* [Distances among centroids](https://learninghub.primer-e.com/link/1046) (in the space of a chosen resemblance measure);
* [Distances from individual points to their own group centroid (in PERMDISP)](https://learninghub.primer-e.com/link/272); and
* [Monte Carlo p-values (in PERMANOVA)](https://learninghub.primer-e.com/link/242).
PCO axes are also used (obviously!) to perform [canonical analysis of principal co-ordinates (CAP)](https://learninghub.primer-e.com/books/permanova-for-primer-guide-to-software-and-statistical-methods/chapter/chapter-5-canonical-analysis-of-principal-coordinates-cap).
---
¶*It is sometimes stated that PCO and metric MDS are the same thing. This is simply incorrect and the example here clearly demonstrates this. PCO is a projection onto principal axes, whereas metric MDS is distance-preserving and minimises stress for a configuration drawn in a chosen number of dimensions. PCO may sometimes be referred to as 'classical scaling', but it should not be confused with multi-dimensional scaling (MDS), metric or otherwise.*
# 16. Means plots
# 16.1 A plot of means with error bars
#### Overview
A plot of the means (or averages) for groups of samples, with error bars capturing some aspect of each group’s variability, is an important tool in univariate data analysis. For example, one may use PERMANOVA to do univariate analysis of variance on a single variable¶, where the [null hypothesis](https://learninghub.primer-e.com/link/1026#bkmrk-and-our-null-hypothe) is focused on the absence of effects defined by individual factors. This is equivalent to a statement that the population means are equivalent across the groups. The alternative hypothesis is that some (one or more) population means differ from some other (one or more) population means across the groups.
The natural visual accompaniment to such an analysis is a plot of the means for the groups of samples corresponding to levels of the factor. It is natural also to show some measure of the variation around those means; i.e., to achieve some way of assessing the expected variability in those mean values in order to assist in the comparison of their values. Error bars constructed around each of the means can accomplish this. Furthermore, if there is more than one factor in the study design, then a plot of the means (and associated error bars) for groups identified by the levels of one factor, calculated and displayed separately across levels of a second (and/or third) factor, yields a highly desirable visualisation of results.
The purpose of the new **Plots** > **Means Plot...** tool in PRIMER 8 with PERMANOVA+ is to make it extremely easy to produce desirable means plots for one or more univariate variables for one-factor or multi-factor study designs. Specifically, using the **Means Plot** dialog in PRIMER 8, you can display group means (as either coloured symbols or bars) for combinations of levels of factors, and you are given a suite of potential methods for drawing suitable error bars on those means. This is a major expansion and improvement over PRIMER 7, where the **Means Plot...** routine only allowed a single factor to be specified, only showed points for the means, and only displayed a $95$ % confidence interval based on the $t$-statistic. See Fig. 15.1 below for a comparison of the old *vs* the new dialog. (Click the image to expand the view).
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/01-fig-01-p7-vs-p8-means-plot.png)
*Fig. 15.1. Comparison of the 'Means plot...' dialog in PRIMER 7 compared to the much more extensive new dialog available in PRIMER 8.*
#### Error bar options
In PRIMER 8, you can choose from the following options to draw error bars around each mean on the plot:
- ± 1 standard error (default);
- ± 1 standard deviation; or
- a confidence interval (choose the percentage) based on:
- the $t$-distribution (useful if the variance is unknown);
- the normal ($z$-)distribution (useful if the variance is known); or
- bootstrap percentiles (for a specified number of bootstrap re-samples), done separately within each group.
For all of the methods apart from the bootstrap percentile method, you can choose to estimate a separate variance for each group, or to estimate a common variance (pooled across the groups).
### Confidence intervals
Consider a random variable $Y$ of unknown distribution (not necessarily normal) having a mean $\mu$ and variance $\sigma^2$. Let $\lbrace y_1, y_2, ..., y_n \rbrace$ be a random sample of $n$ values from this distribution. Let $\bar{y} = \sum_{i=1}^n y_i$ be the mean of the sample. The ***central limit theorem*** states that the distribution of means, $\bar{Y}$, obtained under repeated random sampling in this way, will be normal, with mean $\mu$ and variance $\sigma^2/n$. The standard error, $se$ = $\sqrt{s^2/n}$, where $s^2$ is the [sample variance](https://learninghub.primer-e.com/link/957#bkmrk-average%3A-%24%5Chspace%7B1m), is an unbiased estimator of the ***standard deviation of this distribution of means*** (i.e., $\sigma_{\bar{y}}$). Thus, a plot of $\bar{y} \pm 1 \times se$ is a perfectly natural way to show the variability in the mean calculated from the sample, regardless of that particular variable’s underlying distribution (whether it be normal or not).
If you choose to construct a confidence interval ($\text{CI}$) based on either the $t$-distribution or the normal ($z$)-distribution, then the error bars will be drawn to create a symmetric interval around the mean by multiplying the appropriate quantile from the distribution chosen (e.g., for a $95$ % $\text{CI}$ using the $t$-distribution, the quantile is $1.96$) multiplied by the estimated standard error ($se$) for that group.
Suppose you choose to construct a confidence interval of $\text{CI} = (1-\alpha)\times 100$ % for the mean of a population using the $t$-distribution (e.g., if $\alpha = 0.05$, then the $\text{CI}$ is $95$ %). The interpretation of the interval, under the central limit theorem, is this: under repeated sampling of the population, you would expect $(1-\alpha)\times 100$ % of the confidence intervals constructed in this way to contain the true mean of the population ($\mu$).
### Bootstrap percentiles
The bootstrap percentile error bars are drawn for a given group as follows. Let $n_{b}$ be the requested number of bootstrap samples.
1. Obtain a bootstrap sample by sampling the data values with replacement from the group.
2. Calculate the mean of the bootstrap sample.
3. Repeat steps 1 and 2 a total of $n_b$ times to generate an empirical distribution of bootstrap means.
4. Obtain the following quantiles from the bootstrap distribution: $\alpha/2$ and $1-\alpha/2$, where $\alpha = (100 - \text{CI}) / 100$. Thus, for $\text{CI} = 95$ %, we would get the $0.025$ and $0.975$ quantiles from the bootstrap distribution of means to draw the error bars.
*Note:* If the total number of unique sets of values that can be obtained under bootstrap re-sampling is less than the requested number of bootstraps, then the calculation of the percentile values for the error bars will be done using that total number instead of the requested $n_b$. The formula for the total number of unique sets for a given group is $C(2n-1, n)$, where $C()$ is the binomial coefficient function, and $n$ is the number of samples in the group.
Error bars drawn using bootstrap percentiles have the advantage of being constructed empirically from the distribution of the sample values. They are therefore not necessarily going to be symmetric around the mean. Another thing to note is that although a bootstrap percentile interval constructed in this way *does* give an indication of the variability in the means under repeated sampling, it does *not* have the same interpretation as a classical confidence interval built on the basis of the central limit theorem. The bootstrap is known to have a negative bias in the estimation of the variance of a sample (see, for example, {{@954#bkmrk-andersonetal2017}}). Thus, the bootstrap percentile interval will under-estimate the true distance between the upper and lower percentiles in the underlying population distribution of means that would be obtained under reepeated sampling. No bias-correction has been implemented in the construction of the bootstrap percentile error bars in PRIMER 8 for means plots; however, if the sample size ($n$) is reasonably large, then the bias will be relatively small (of order $1/n$)†.
---
¶ *This is achieved by basing the PERMANOVA analysis on a Euclidean distance matrix calculated from a single variable.*
†*In fact, in Appendix B of {{@954#bkmrk-andersonetal2017}}, it is shown that the exact downward bias in the bootstrap estimate of the variance of the mean of a univariate random variable is $(1−1/n)$.*
# 16.2 Example: Fal biota (one-way case)
Let's look at an example of a means plot for a one-way case using biotic data ($p$ = 131 taxa) from the Fal estuary. The data consist of counts of organisms living in benthic sediments obtained from $n$ = 5-7 sites in each of 5 creeks running into the Fal estuary, Cornwall, UK ({{@954#bkmrk-somerfieldetal1994a}}, {{@954#bkmrk-somerfieldetal1994b}}). These 5 creeks had differing levels of heavy metals in their sediments, due to historical tin and copper mining in their respective valleys. Interest lies in identifying whether the 5 creeks also differed in their diversity of benthic infauna. The full set of biotic data (macrofauna and also the meiofauna - comprised of nematodes and copepods) are contained in the file 'Fal_all_biotic_taxa.pri', found in the 'Fal_benthic_fauna' folder in 'Examples_P8'.
1. Start by calculating some univariate diversity measures for each sample. Open the 'Fal_all_biotic_taxa.pri' file in PRIMER 8, then click **Analyse** > **DIVERSE...**, take all the defaults and click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/02-diverse-fal-estuary.png)
This will produce a sheet, called 'Data1' with values for various diversity indices for each of the sampling units.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/03-diversity-data-fal.png)
2. Calculate a means plot that shows the means $\pm$ 1 standard error for the 5 creeks for each of these diversity measures. From the 'Data1' sheet, click **Plots** > **Means Plot...**.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/04-diverse-fal-means-plot-menu.png)
Leave all of the defaults and click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/05-diverse-fal-means-plot-dialog-b.png)
The result will be a 'Multi plot' item, that includes 6 individual plots - one for each of the 6 diversity measure variables in the original 'Data1' sheet.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/06-diverse-fal-multi-plot.png)
3. View an individual plot by clicking on it. For example, you can click on the plot in the top left-hand corner of 'MultiPlot1' to see a means plot for the variable of species richness ('S', the number of different taxa). Alternatively, you can click on 'Graph1' in the Explorer tree.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/07-diverse-fal-species-richness.png)
4. To change the colours, symbols or line types, you can click **Graph** > **Special...**, which shows you a dialog with a plethora of options to modify all of these details of your means plot. For example, you can click the 'Key...' button next to the factor to change the colours/symbols, or choose to $\bullet$ **Join means**, including the option to choose a custom colour for the joining line, like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/08-diverse-fal-s-join-means.png)
This will yield the following modified plot:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/09-diverse-fal-s-plot-with-line.png)
5. You can also choose to show the means as bars instead. Click **Graph** > **Special...**, then (Display Means as > $\bullet$ Bars) and (Error bars > Appearance > $\bullet$ Custom).
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/10a-diverse-fal-s-barscustom-dialog.png)
This yields the following plot:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/10-diverse-fal-s-bars.png)
Next, we shall consider some multi-factor plots.
# 16.3 Example: Okura macrofauna (two-way nested case)
Let's look at a means plot where we have a two-way nested design. We will examine data consisting of counts of benthic macrofauna from intertidal sites in the Okura estuary, located north of Auckland, New Zealand ({{@954#bkmrk-andersonetal2004}}). Hydrodynamic models were used to identify different areas within the estuary as having high (H), medium (M) or low (L) probabilities of sediment deposition. There were several areas of each type (H, M, or L) interspersed along the estuary from its mouth to its inner reaches. Interest lies in identifying potential differences in the abundances of infaunal species at sites having different probabilities of sediment deposition. Some organisms are expected to be more tolerant of sediment inputs, while others (particularly filter feeders) may be more vulnerable to fine sediment inputs.
Six (6) sediment cores (13 cm in diameter × 15 cm deep) were obtained from random positions within each of 15 sites along the estuary, with 5 sites from each of the high, medium and low depositional types of environments. Sampling was repeated 6 times (twice in each of three seasons in 2001-2002), yielding a total of 36 cores sampled per site.
These data are contained in the file 'Okura_macrofauna.pri', found in the 'Okura_macrofauna' folder in 'Examples_P8'. We shall ignore, for now, the temporal factors in the study and focus only on the spatial factors in a two-way nested design:
- Deposition (fixed with $a$ = 3 levels: H, M or L); and
- Site (random and nested, $b$ = 5 sites within each Deposition, thus 15 sites in total)
We will create a means plot that shows the mean ± 1 standard error for each site for a single univariate variable: the abundance of cockles, *Austrovenus stutchburyi*.
1. Open the 'Okura_macrofauna.pri' file in PRIMER 8, look for the variable (column) named 'Austrovenus stutchburyi', highlight this species by clicking on its name, then click **Select** > **Highlighted**.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/11-okura-select-austrovenus.png)
2. You will see the selected species as a single column with a blue background. With this single species selected, click **Plots** > **Means Plot...**, like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/12-okura-means-plot-menu.png)
3. Given that the design is nested, we would like to calculate a separate mean (and standard error) for each site, and we would like bars corresponding to different sites to be arranged along the x-axis separately within each of the deposition groups (H, M and L). We therefore need to choose the following in the 'Means Plot' dialog:
- Factor A (different symbols/colours) > Site
- $\checkmark$ Split into separate groups (along the x-axis) by Factor B > Deposition
- Display means as > $\bullet$ Bars
as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/13-okura-means-plot-dialog.png)
Click '**OK**'.
The resulting plot (given in the output file called 'Graph1' initially will simply show a different coloured bar for each of the sites, like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/14-okura-initial-plot.png)
Although this is not too bad, the plot could be improved in a number of ways. First, notice that the 3 levels of 'Deposition' are ordered alphabetically, but we might prefer them to have a more logical order, e.g. L < M < H. Second, we might like to make the error bars black, rather than being the same colour as the bars themselves, so that we can more easily see their extent above and below the mean value. Third, it might be nice if the colours of the bars corresponded to the three different levels of 'Deposition', rather than being a different colour for every site. After all, 'Site' is a random factor.
#### Special graphical properties of Means Plots
Means plots in PRIMER have a number of special graphical properties that can be altered by you even after the plot has been made. Click **Graph** > **Special...** to see them. Specifically, you can choose:
- to change the order of factor levels or the colours/symbols/lines associated with particular factors (*via* the 'Key' menu)
- to display means as points or bars
- whether and how the means may be joined, and the associated line types for joining
- the sizes of gaps within or between groups of means
- whether to show separator lines between groups and what line colour/type to use for this
- the colour and style of error bar lines, and whether to show both sides or just one side (the outer side).
4. To change the colour of the error bars in the plot of *Austrovenus* ('Graph1'), click **Graph** > **Special...** and under 'Error bars', choose 'Appearance > $\bullet$ Custom', and choose a colour of your choice (the default of black will do fine here). Next, to change the order of the deposition groups on the X axis, click the '**Key**' button next to 'Factor B: Deposition'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/15-means-plot-special-dialog.png)
Inside the 'Key' menu, click on the level labeled 'H', and move it progressively to the bottom using the 'Move' arrow button [](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/image-1784168270813.png), so that the order of the Deposition levels (from top to bottom) is changed to 'L', followed by 'M', followed by 'H', respectively, then click '**OK**', as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/15c-key-final.png)
This yields a plot that is better, viz:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/16-okura-modified-plot.png).
5. We have improved the colour of the error bars and the ordering of the factor levels along the x axis. Now let's change the colours of the bars to match the factor of Deposition. Click **Graph** > **Sample Labels & Symbols...** and choose the
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/17-graph-options-dialog.png)
The final plot looks like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/18-okura-final-plot.png)
Here, we can clearly see both the substantial variability among sites in the mean abundance of cockles. We can also see that, on average, the mean abundance of cockles is greater at sites having either a low or medium probability of sediment deposition (amber or green bars), compared to sites having a high probability of sediment deposition (blue bars).
# 16.4 Example: Leschenault fish (two-way crossed case)
An example of a two-way crossed design is provided by a study of trawl samples for fish communities (using a 21.5m seine net) from 4 regions - basal (B), lower (L), upper (U), and apex (A) - occurring between the mouth and upper reaches of the Leschenault estuary, Western Australia ({{@954#bkmrk-vealeetal2014}}). Samples were taken from each region over 4 seasons - spring (Sp), summer (S), autumn (A), and winter (W), with 6-8 replicate samples from each of these 16 combinations (taken over the years 2008-2010). Data are located in the file 'Leschenault_fish_counts.pri', found in the 'Leschenault_fish' folder in 'Examples_P8'. This is a simplified version of the data described by {{@954#bkmrk-vealeetal2014}}. The following two factors are crossed with one another, as all regions were sampled in every season:
- **Region**, fixed with 4 levels: basal (B), lower (L), upper (U), and apex (A)
- **Season**, fixed wtih 4 levels: spring (Sp), summer (S), autumn (A), and winter (W)
#### Compare seasonal means for different regions
We shall focus our attention here on just a few species that have high frequencies of occurrence across the dataset as a whole. We wish to construct means plots to examine and compare potential changes in the average abundance of a given fish species across the seasons (considered separately for each region) and across the regions (considered separately for each season).
1. Begin by selecting the top 4 species, based on their percentage contribution to any one sample. Open the file 'Leschenault_fish_counts.pri', found in the 'Leschenault_fish' folder in 'Examples_P8' in PRIMER. Click **Select** > **Variables...** and choose:
$\bullet$ Top 4 variable(s) based on: > $\bullet$ Percent contribution to any one sample, as shown below.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/19-lesch-data-select-top4.png)
2. First we'll make plots that permit easy comparison of the seasonal means, separately within each region. Click **Plots** > **Means Plot...**, then choose the following:
- Factor A (different symbols/colours): Season
- $\checkmark$Split into separate groups (along the x-axis) by Factor B > Region > $\checkmark$Draw group separator lines
- Display means as > $\bullet$ Points
- Join means > $\bullet$ Across levels of Factor A (within levels of Factor B)
Leave the rest as defaults and click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/20-lesch-means-plot-season-dialog.png)
This will produce 4 graphics - one for each species - presented in the multiplot object 'MultiPlot1':
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/21-leschenault-multiplot.png)
3. Let's take a closer look at the fish species '246010' (actual species names were not provided with this particular dataset). Click on 'Graph2' in the Explorer tree (or click on the upper right-hand plot in the multiplot graphic).
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/22-lesch-species-246010.png)
This species generally had greater average abundance in autumn (green symbols), particularly in the upper (U) and apex (A) regions of the estuary.
#### Compare regional means for different seasons
4. We can change the colour or line-type of the joining lines. Note that the default colours and line-types for these are based on the colour key for Factor B, which is 'Region' here. We can also change the joining lines so that they occur across the regions, rather than across the seasons. A further option is to change the colour or style (or even remove entirely) the dotted lines separating the different regions along the x-axis. Let's check out the visual effect of these changes on our graphic. Click **Graph** > **Special...**. Then, in the 'Means Plot' special dialog, choose the following:
- Join means > $\bullet$ Across levels of Factor B (within levels of Factor A)
- Separator lines > $\bullet$ No separate lines
then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/23-lesch-change-joining-lines.png)
The resulting plot will look like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/23-lesch-first-revised-plot.png)
5. Another tweak we might like to consider here is to change the relative gap width among the categories (bars or points) that are being plotted within *vs* between the groups along the x-axis. For example, we could make the seasonal means within a group closer together, with a bigger gap between the different regions. Click **Graph** > **Special...**. Then, in the 'Means Plot' special dialog, in the section entitled 'Choose gaps', make the 'Across group gap' six times (say) the 'Within group gap', e.g.:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/24-choose-gap-width.png)
The resulting gap-adjusted plot is quite a bit clearer, *viz.*:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/25-lesch-2nd-revised-plot.png)
In this graphic, we can now more clearly see there is a gradual increase in the mean abundance of this fish species as you go up the estuary (from the basal to the lower to the upper and finally to the apex region), but this trend is only really apparent in the autumn season (green symbols), and not so apparent in the other seasons.
6. Another possibility would be to ***swap the roles*** of the factors in this crossed design for the plot itself to begin with. We can examine the regional means separately within each season. To do this, go back to the original data sheet (called 'Leschenault_fish_counts' in the Explorer tree) and click **Plots** > **Means Plot...**, then choose 'Region' as Factor A and 'Season' as Factor B, leaving all the rest as before.
A Multiplot will be produced (4 graphics for the 4 species), with the top-right plot ('Graph6') corresponding once again to fish species 246010, like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/26-lesch-3rd-revised-plot.png)
The general trend of increasing mean abundance as you go from the basal to the apex region of the estuary (as previously described) is even more plainly seen here, and is very marked in autumn.
Note also that the joining lines here correspond to the underlying colours/line-types of the regions. If we wish, we can make these joining lines a solid colour of our choice. For example, we shall make them a solid charcoal colour. Click **Graph** > **Special...**, then under 'Line appearance' in the dialog, choose '$\bullet$ Custom' and click on the solid box next to the word 'Colour:' and you can choose your colour explicitly, like so:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/27-line-colour-all.png)
The resulting plot that uses a solid charcoal colour for the joining lines looks like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/28-lesch-fourth-revised-plot.png)
#### Bars instead of points
8. Of course, we can also look at this plot using bars instead of points. Simply click **Graph** > **Special...**, then choose:
- Display means as > $\bullet$ Bars
- Error bars > Appearance > $\bullet$ Custom
to get the following plot:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/29-lesch-as-bar-plot.png)
It is clear there are many ways to customise and refine these means plots. For example, the colours associated with the individual factor levels can be changed using the 'Key' button(s) in the **Graph** > **Special...** menu dialog. In addition, a great host of other aspects of the plot can also obviously be customised using the **Graph** > **General...** menu dialog, including the text and fonts for titles, the legend ('Keys'), the 'History' text pane, the X or Y axes, and so on.
# 16.5 Example: New Zealand snapper (three-way case)
The **Means Plot** procedure in PRIMER 8 will operate just as other univariate plots do (such as histograms and dotplots) in that multiple plots will be produced - one for each variable - if there are multiple variables in the original spreadsheet. However, means plots also provide the end-user with the ability to create a multi-plot object by reference to an additional factor. This can be a second or a third factor, depending on the situation. Here, we will show a couple of examples that use this approach advantageously to show patterns across all combinations of three factors in a single multi-plot graphic.
The example dataset here are counts of snapper (*Chrysophrys auratus*) taken from baited remote underwater video (BRUV) apparatuses deployed at each of three locations: Leigh, Tawharanui and Hahei, in each of two seasons over several years (1997-2007), from sites that were located either inside or outside of a marine reserve (the factor here is 'Status', with levels 'NR' = non-reserve and 'R' = reserve). See {{@954#bkmrk-smithetal2014}} for details.
#### Status by Year by Location
First, we will visualise patterns in mean relative abundances (total count per BRUV unit) of snapper inside *vs* outside reserves across the years at the different locations. The data are located in the file named 'NZ_Snapper_counts.pri', located in the folder 'Examples_P8' > NZ_Snapper_counts'.
1. Open the file 'NZ_Snapper_counts.pri' in PRIMER and click on the column (variable) name 'tot.snapper' in the sheet to highlight it, then click **Select** > **Highlighted**. The selected variable will be highlighted in blue, thusly:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/30-snapper-counts.png)
Next, we shall plot the mean count of total snapper observed (per BRUV) with along with the 2.5% and 97.5% percentiles of the distribution of bootstrap means from 1000 bootstrap re-samples.
2. Click **Plots** > **Means Plot...** and choose the following options in the 'Means plot' dialog:
- Factor A (different symbols/colours) > Status
- $\checkmark$Split into separate groups (along the x-axis) by Factor B > Year, and 'untick' the box '$\Box$ Draw group separator lines'
- $\checkmark$Split into separate plots by Factor C > Location
- Display means as > $\bullet$ Points
- Join means > $\bullet$ Across levels of Factor B (within levels of Factor A)
- Error bars > $\bullet$ Bootstrap percentile
as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/31-snapper-means-plot-dialog-1.png)
The resulting multi-plot shows the means (and bootstrap percentiles) for snapper counts from reserve and non-reserve sites separately in each year, with different locations being shown as different graphics within the multi-plot object, *viz*:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/32-snapper-multiplot.png)
Examining the plot for Leigh only (click on the upper right-hand plot, or click 'Graph2' in the Explorer tree), we see, for example:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/33-snapper-leigh-only.png)
You can change or arrange the graphics in a multi-plot (i.e., by specifying the number of columns and rows of graphics that you want, and also their ordering) etc. by first clicking on the multiplot item itself in the Explorer tree (e.g., 'MultiPlot1'), then click **Graph** > **Special...**.
#### Season by Status by Location
We might also like to look at differences in the mean relative abundances of snapper inside *vs* outside these reserves, also split by location. We will do this just for Leigh and Hahei, where we have more years of information.
4. From the 'NZ_Snapper_counts' data sheet, click **Select** > **Samples...**, then choose $\bullet$ Factor levels > Location and click the 'Levels...' button [](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/levels-button.png). Move the words 'Hahei' and 'Leigh' from the 'Available' box on the left to the 'Include' box on the right (by clicking on them each and then clicking on the right-arrow button [](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/right-arrow.png)), then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/35-select-locations.png)
5. Click **Plots** > **Means Plot...** and choose the following options in the 'Means plot' dialog:
- Factor A (different symbols/colours) > Season
- $\checkmark$Split into separate groups (along the x-axis) by Factor B > Status, and $\checkmark$ Draw group separator lines
- $\checkmark$Split into separate plots by Factor C > Location
- Display means as > $\bullet$ Bars
- Join means > $\bullet$ Across levels of Factor B (within levels of Factor A)
- Error bars > $\bullet$ Bootstrap percentile
as shown below:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/34-snapper-seasons-by-status.png)
The resulting graphic, by default, looks like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/36-snapper-season-graphic.png)
6. We might like to change the colours of these bars to something a bit more reflective of spring and autumn. Go back to the original data sheetcalled 'NZ_Snapper_counts' and click **Edit** > **Factors...**, then click on the column called 'Season', and click the 'Key' button [](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/key-button.png). In the 'Key' dialog, click on the colour box for 'Spring' and choose (say) a green colour, then click on the colour box for 'Autumn' and choose (say) an autumnal colour of brown or orange/red, then click '**OK**'.
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/37-edit-key.png)
Re-running the same dialog for the means plot now as we ran before (see step 5 above), we obtain the following ('MultiPlot3'):
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/38-snapper-season-graphic-revised.png)
7. We might also like to change the error bar colours so that they are black, making it easier to see the lower and upper bootstrap percentile values. For each of 'Graph6' and 'Graph7', in turn (both of these graphics are in 'MultiPlot3'), click **Graph** > **Special...**, then choose
Error bars > Appearance > $\bullet$ Custom (and leave the default colour of black here).
The resulting plot looks like this:
[](https://learninghub.primer-e.com/uploads/images/gallery/2026-07/39-snapper-season-graphic-revised2.png)
The above plot shows that the mean relative abundance of snapper (sampled using BRUVs) is typically greater in the autumn than in the spring, on average. Also, the largest seasonal difference in mean was observed inside the reserve at Leigh.
#### Means plots to accompany PERMANOVA formal analysis
Typcially, we would perform a more complete univariate statistical analysis of these data using (say) PERMANOVA in accordance with the full multi-factor sampling design. This could be done easily in PRIMER with PERMANOVA+ by calculating Euclidean distances among the samples based on the single variable of 'tot.snapper' and creating an appropriate Design file corresponding to the full multi-factor experimental design (including the nested factor of 'Areas' as well). Pair-wise comparisons could then be run to follow up any significant interactions that might be discovered among the main factors of interest.
It is useful to show relevant multi-factor means plots to accompany PERMANOVA tests and pair-wise comparisons. For example, a statistically significant threee-way interaction of Season$\times$Status$\times$Location might well be accompanied by not only the relevant pair-wise comparisons, but also the above means plots that permit one to compare the means (and their variability) visually in a useful way.
#### Output means data to worksheet
There are many ways to customise means plots within PRIMER. It is also possible (in the 'Means Plot' dialog) to tick the box to '$\checkmark$Output means data to worksheet'. This will output the means and standard errors (or the means and the upper and lower bounds of error bars of your choice) which are calculated to produce the plot. If there are customisations that you want to create in a different graphics package, or if you would rather get the values for the means (and their measures of variation) in a tabled format (e.g., for publication in supplementary material, etc.), then this extra tool is very handy, and may well be more efficient than using **Tools** > **Summary Stats...**, even though the latter has been [expanded considerably](https://learninghub.primer-e.com/books/whats-new-in-primer-8/chapter/1-expanded-summary-statistics) from what was available in version 7.
# References
---
[ Adegoke (2019) ](https://learninghub.primer-e.com/link/954#bkmrk-adegoke2019)
Adegoke, N. (2019) *Contributions to improve power, efficiency and scope of control-chart methods* . PhD thesis, School of Natural and Computational Sciences (SNCS), Massey University, New Zealand, 284 pp.
[link to source](https://mro.massey.ac.nz/server/api/core/bitstreams/5ae914c9-bdc1-4d11-982d-fc299245e2e0/content)
---
[ Adegoke *et al.* (2018) ](https://learninghub.primer-e.com/link/954#bkmrk-adegokeetal2018)
Adegoke, N.A., Smith, A.N.H., Anderson, M.J., Abbasi, S.A. & Pawley, M.D.M. (2012) Shrinkage estimates of covariance matrices to improve the performance of multivariate cumulative sum control charts. *Computers & Industrial Engineering*, **117**, 207-216.
---
[ Ahmad (2014) ](https://learninghub.primer-e.com/link/954#bkmrk-ahmad2014)
Ahmad, M.R. (2014) A *U*-statistic approach for a high-dimensional two-sample mean testing problem under non-normality and Behrens–Fisher setting. *Annals of the Institute of Statistical Mathematics*, **66**, 33-61.
---
[ Ahmad *et al.* (2012) ](https://learninghub.primer-e.com/link/954#bkmrk-ahmadetal2012)
Ahmad, M.R., vonRosen, D. & Singull, M. (2012) A note on mean testing for high dimensional multivariate data under non-normality. *Statistica Neerlandica*, **67**, 88-99.
---
[ Anderson (1992) ](https://learninghub.primer-e.com/link/954#bkmrk-anderson1992)
Anderson, M.J. (1992) *Settlement and succession of oysters and fouling organisms in Quibray Bay, New South Wales*. Graduate Diploma of Science (Honours) thesis, School of Biological Sciences, University of Sydney, NSW, Australia, 51 pp.
---
[ Anderson (1996) ](https://learninghub.primer-e.com/link/954#bkmrk-anderson1996)
Anderson, M.J. (1996) A chemical cue induces settlement of Sydney Rock oysters, *Saccostrea commercialis*, in the laboratory and in the field. *Biological Bulletin*, **190**, 350-358.
---
[ Anderson (2001) ](https://learninghub.primer-e.com/link/954#bkmrk-anderson2001)
Anderson, M.J. (2001) A new method for non-parametric multivariate analysis of variance. *Austral Ecology*, **26**, 32-46.
---
[ Anderson (2006) ](https://learninghub.primer-e.com/link/954#bkmrk-anderson2006)
Anderson, M.J. (2006) Distance-based tests for homogeneity of multivariate dispersions. *Biometrics*, **62**, 245-253.
---
[ Anderson (2017) ](https://learninghub.primer-e.com/link/954#bkmrk-anderson2017)
Anderson, M.J. (2017) Permutational multivariate analysis of variance (PERMANOVA). *Wiley StatsRef: Statistics Reference Online*, **stat07841**, 15 pp.
[link to source](https://onlinelibrary.wiley.com/doi/full/10.1002/9781118445112.stat07841)
---
[ Anderson *et al.* (2005) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2005)
Anderson, M.J., Diebel, C.E., Blom, W.M. & Landers, T.J. (2005) Consistency and variation in kelp holdfast assemblages: spatial patterns of biodiversity for the major phyla at different taxonomic resolutions. *Journal of Experimental Marine Biology and Ecology*, **320**, 35-56.
---
[ Anderson *et al.* (2006) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2006)
Anderson, M.J., Ellingsen, K.E. & McArdle, B.H. (2006) Multivariate dispersion as a measure of beta diversity. *Ecology Letters*, **9**, 683-693.
---
[ Anderson *et al.* (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2004)
Anderson, M.J., Ford, R.B., Feary, D.A. & Honeywill, C. (2004) Quantitative measures of sedimentation in an estuarine system and its relationship with intertidal soft-sediment infauna. *Marine Ecology Progress Series*, **272**, 33-48.
---
[ Anderson *et al.* (2025) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2025)
Anderson, M.J., Gorley, R.N. & Terlizzi, A. (2025) The incremental progression from fixed to random factors in the analysis of variance: a new synthesis. *Australian & New Zealand Journal of Statistics*, **67**, 3-30.
---
[ Anderson & Legendre (1999) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonlegendre1999)
Anderson, M.J. & Legendre, P. (1999) An empirical comparison of permutation methods for tests of partial regression coefficients in a linear model. *Journal of Statistical Computation and Simulation*, **62**, 271-303.
---
[ Anderson & Millar (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonmillar2004)
Anderson, M.J. & Millar, R.B. (2004) Spatial variation and effects of habitat on temperate reef fish assemblages in northeastern New Zealand. *Journal of Experimental Marine Biology and Ecology*, **305**, 191-221.
---
[ Anderson & Robinson (2001) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonrobinson2001)
Anderson, M.J. & Robinson, J. (2001) Permutation tests for linear models. *Australian & New Zealand Journal of Statistics*, **43**, 75-88.
---
[ Anderson & Thompson (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonthompson2004)
Anderson, M.J. & Thompson, A.A. (2004) Multivariate control charts for ecological and environmental monitoring. *Ecological Applications*, **14**, 1921-1935.
---
[ Anderson *et al.* (2013) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2013)
Anderson, M.J., Tolimieri, N. & Millar, R.B. (2013) Beta diversity of demersal fish assemblages in the north-eastern Pacific: interactions of latitude and depth. *PLoS ONE*, **8(3)**, e57918.
[link to source](https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0057918)
---
[ Anderson & Underwood (1994) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonunderwood1994)
Anderson, M.J. & Underwood, A.J. (1994) Effects of substratum on the recruitment and development of an intertidal estuarine fouling assemblage. *Journal of Experimental Marine Biology and Ecology*, **184**, 217-236.
---
[ Anderson & Walsh (2013) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonwalsh2013)
Anderson, M.J. & Walsh, D.C.I. (2013) PERMANOVA, ANOSIM, and the Mantel test in the face of heterogeneous dispersions: what null hypothesis are you testing? *Ecological Monographs*, **83**, 557-574.
---
[ Anderson *et al.* (2017) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2017)
Anderson, M.J., Walsh, D.C.I., Clarke, K.R., Gorley, R.N. & Guerra-Castro, E. (2017) Some solutions to the multivariate Behrens-Fisher problem for dissimilarity-based analyses. *Australian & New Zealand Journal of Statistics*, **59**, 57-79.
---
[ Anderson *et al.* (2022) ](https://learninghub.primer-e.com/link/954#bkmrk-andersonetal2022)
Anderson, M.J., Walsh, D.C.I., Sweatman, W.L. & Punnett, A.J. (2022) Non-linear models of species' responses to environmental and spatial gradients. *Ecology Letters*, **25**, 2739-2752.
[link to source](https://onlinelibrary.wiley.com/doi/10.1111/ele.14121)
---
[ Barnard (1959) ](https://learninghub.primer-e.com/link/954#bkmrk-barnard1959)
Barnard, G.A. (1959) Control charts and stochastic processes. *Journal of the Royal Statistical Society, Series B*, **21**, 239-271.
---
[ Behrens (1929) ](https://learninghub.primer-e.com/link/954#bkmrk-behrens1929)
Behrens, W.V. (1929) Ein Beitrag zur Fehlerberechnung bei wenigen Beobachtungen. *Landwirtschaftliches Jaresbch*, **68**, 807–837.
---
[ Belloni & Didier (2008) ](https://learninghub.primer-e.com/link/954#bkmrk-bellonididier2008)
Belloni, A. & Didier, G. (2008) On the Behrens-Fisher problem: a globally convergent algorithm and a finite sample study of the Wald, LR and LM tests. *The Annals of Statistics*, **36**, 2377-2408.
---
[ Bliss (1958) ](https://learninghub.primer-e.com/link/954#bkmrk-bliss1958)
Bliss, C.I. (1958) Periodic regression in biology and climatology. *The Connecticut Agricultural Experiment Station*, Bulletin **615**, 1-55.
---
[ Boik (1987) ](https://learninghub.primer-e.com/link/954#bkmrk-boik1987)
Boik, R.J. (1987) The Fisher-Pitman permutation test: anon-robut alternative to the normal theory F test when variances are heterogeneous. *British Journal of Mathematical and Statistical Psychology*, **40**, 26-42.
---
[ Borg & Groenen (2005) ](https://learninghub.primer-e.com/link/954#bkmrk-borggroenen2005)
Borg, I. & Groenen, P.J.F. (2005) *Modern multidimensional scaling, 2nd edition*. New York, NY, USA: Springer.
---
[ Box (1954) ](https://learninghub.primer-e.com/link/954#bkmrk-box1954)
Box, G.E.P. (1954) Some theorems on quadratic forms applied in the study of analysis of variance problems, I: effect of inequality of variance in the one-way classification. *The Annals of Mathematical Statistics*, **25**, 290-302.
---
[ Brown & Forsythe (1974) ](https://learninghub.primer-e.com/link/954#bkmrk-brownforsythe1974)
Brown, M.B. & Forsythe, A.B. (1974) The small sample behaviour of some statistics which test the equality of several means. *Technometrics*, **16**, 129-132.
---
[ Christensen & Rencher (1997) ](https://learninghub.primer-e.com/link/954#bkmrk-christensenrencher1997)
Christensen, W.F. & Rencher, A.C. (1997) A comparison of type I error rates and power levels for seven solutions to the multivariate Behrens-Fisher problem. *Communications in Statistics - Simulation and Computation*, **26**, 1251-1273.
---
[ Clarke (1993) ](https://learninghub.primer-e.com/link/954#bkmrk-clarke1993)
Clarke, K.R. (1993) Nonparametric multivariate analyses of changes in community structure. *Australian Journal of Ecology*, **18**, 117-143.
---
[ Clarke & Ainsworth (1993) ](https://learninghub.primer-e.com/link/954#bkmrk-clarkeainsworth1993)
Clarke, K.R. & Ainsworth, M. (1993) A method of linking multivariate community structure to environmental variables. *Marine Ecology Progress Series*, **92**, 205-219.
---
[ Clarke *et al.* (2006a) ](https://learninghub.primer-e.com/link/954#bkmrk-clarkeetal2006a)
Clarke, K.R., Chapman, M.G., Somerfield, P.J. & Needham, H.R. (2006a) Dispersion-based weighting of species counts in assemblage analyses. *Marine Ecology Progress Series*, **320**, 11-27.
---
[ Clarke *et al.* (2006b) ](https://learninghub.primer-e.com/link/954#bkmrk-clarkeetal2006b)
Clarke, K.R., Somerfield, P.J. & Chapman, M.G. (2006b) On resemblance measures for ecological studies, including taxonomic dissimilarities and a zero-adjusted Bray-Curtis coefficient for denuded assemblages. *Journal of Experimental Marine Biology and Ecology*, **330**, 55-80.
---
[ Clarke *et al.* (2014) ](https://learninghub.primer-e.com/link/954#bkmrk-clarkeetal2014)
Clarke, K.R., Gorley, R.N., Somerfield, P.J. & Warwick, R.M. (2014) *Change in marine communities: an approach to statistical analysis and interpretation, 3rd edition*. Plymouth, UK: PRIMER-E.
[link to source](https://learninghub.primer-e.com/books/change-in-marine-communities)
---
[ Clinch & Keselman (1982) ](https://learninghub.primer-e.com/link/954#bkmrk-clinchkeselman1982)
Clinch, J.J. & Keselman, H.T. (1982) Parametric alternatives to the analysis of variance. *Journal of Educational Statistics*, **7**, 207-214.
---
[ Conover (1972) ](https://learninghub.primer-e.com/link/954#bkmrk-conover1972)
Conover, W.J. (1972) On methods of handling ties in the Wilcoxon signed-rank test. *Journal of the American Statistical Association*, **68**, 985-988.
---
[ Coombs & Algina (1996) ](https://learninghub.primer-e.com/link/954#bkmrk-coombsalgina1996)
Coombs, W.T. & Algina, J. (1996) New test statistics for MANOVA/descriptive discriminant analysis. *Educational and Psychological Measurement*, **56**, 382-402.
---
[ Cornfield & Tukey (1956) ](https://learninghub.primer-e.com/link/954#bkmrk-cornfieldtukey1956)
Cornfield, J. & Tukey, J.W. (1956) Average values of mean squares in factorials. *The Annals of Mathematical Statistics*, **27**, 907-949.
---
[ Danielidis (1991) ](https://learninghub.primer-e.com/link/954#bkmrk-danielidis1991)
Danielidis, D.B. (1991) *A systematic and ecological study of diatoms of the lagoons of Messolongi, Aitoliko and Kleissova (Greece).* PhD thesis, University of Athens, Greece.
---
[ Darling (1957) ](https://learninghub.primer-e.com/link/954#bkmrk-darling1957)
Darling, D.A. (1957) The Kolmogorov-Smirnov, Crámer-von Mises tests. *The Annals of Mathematical Statistics*, **28**, 823-838.
---
[ Dowling (2021) ](https://learninghub.primer-e.com/link/954#bkmrk-dowling2021)
Dowling, A. (2021) *Environmental influence on size frequency distributions of the Pacific blue mussel (Mytilus trossulus) in two glacially influenced estuaries.* MSc thesis, University of Alaska, Fairbanks, AK, USA.
---
[ Ellingsen & Gray (2002) ](https://learninghub.primer-e.com/link/954#bkmrk-ellingsengray2002)
Ellingsen, K.E. & Gray, J.S. (2002) Spatial patterns of benthic diversity: is there a latitudinal gradient along the Norwegian continental shelf? *Journal of Animal Ecology*, **71**, 373-389.
---
[ Efron & Tibshirani (1993) ](https://learninghub.primer-e.com/link/954#bkmrk-efrontibshirani1993)
Efron, B. & Tibshirani, R.J. (1993) *An introduction to the bootstrap*. New York, NY, USA: Chapman and Hall.
---
[ Fisher (1935) ](https://learninghub.primer-e.com/link/954#bkmrk-fisher1935)
Fisher, R.A. (1935) The fiducial argument in statistical inference. *Annals of Eugenics*, **6**, 391-398.
---
[ Friedman (1989) ](https://learninghub.primer-e.com/link/954#bkmrk-friedman1989)
Friedman, J. (1989) Regularized discriminant analysis. *Journal of the American Statistical Association*, **84**, 165-175.
---
[ Gamage *et al.* (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-gamageetal2004)
Gamage, J., Mathew, T. & Weerhandi, S. (2004) Generalized p-values and generalized confidence regions for the multivariate Behrens-Fisher problem and MANOVA. *Journal of Multivariate Analysis*, **8**, 177-189.
---
[ Ghosh & Kim (2001) ](https://learninghub.primer-e.com/link/954#bkmrk-ghoshkim2001)
Ghosh, M. & Kim, Y.-Y. (2001) The Behrens-Fisher problem revisited: a Bayes-frequentist synthesis. *Canadian Journal of Statistics*, **29**, 5-17.
---
[ Glasby (1997) ](https://learninghub.primer-e.com/link/954#bkmrk-glasby1997)
Glasby, T.M. (1997) Analysing data from post-impact studies using asymmetrical analyses of variance: a case study of epibiota on marinas. *Australian Journal of Ecology*, **22**, 448-459.
---
[ Glasby (1999) ](https://learninghub.primer-e.com/link/954#bkmrk-glasby1999)
Glasby, T.M. (1999) Interactive effects of shading and proximity to the seafloor on the development of subtidal epibiotic assemblages. *Marine Ecology Progress Series*, **190**, 113-124.
---
[ Glasby & Underwood (1998) ](https://learninghub.primer-e.com/link/954#bkmrk-glasbyunderwood1998)
Glasby, T.M. & Underwood, A.J. (1998) Determining positions for control locations in environmental studies of estuarine marinas. *Marine Ecology Progress Series*, **171**, 1-14.
---
[ Gower (1966) ](https://learninghub.primer-e.com/link/954#bkmrk-gower1966)
Gower, J.C. (1966) Some distance properties of latent root and vector methods used in multivariate analysis. *Biometrika*, **53**, 325-338.
---
[ Glass *et al.* (1972) ](https://learninghub.primer-e.com/link/954#bkmrk-glassetal1972)
Glass, G.V., Peckham, P.D. & Sanders, J.R. (1972) Consequences of failure to meet assumptions underlying fixed effects analyses of variance and covariance. *Review of Educational Research*, **42**, 237-288.
---
[ Gray *et al.* (1990) ](https://learninghub.primer-e.com/link/954#bkmrk-grayetal1990)
Gray, J.S., Clarke, K.R., Warwick, R.M. & Hobbs, G. (1990) Detection of initial effects of pollution on marine benthos: an example from the Ekofisk and Eldfisk oilfields, North Sea. *Marine Ecology Progress Series*, **66**, 285-299.
---
[ Hartley (1967) ](https://learninghub.primer-e.com/link/954#bkmrk-hartley1967)
Hartley, H.O. (1967) Expectations, variances and covariances of ANOVA mean squares by 'synthesis'. *Biometrics*, **23**, 105-114.
---
[ Hartley *et al.* (1978) ](https://learninghub.primer-e.com/link/954#bkmrk-hartleyetal1978)
Hartley, H.O., Rao, J.N.K. & Lamotte, L.R. (1978) A simple 'synthesis'-based method of variance component estimation. *Biometrics*, **34**, 233-242.
---
[ Hayes (1996) ](https://learninghub.primer-e.com/link/954#bkmrk-hayes1996)
Hayes, A.F. (1997) Permutation test is not distribution free. *Psychological Methods*, **1**, 184-198.
---
[ Horsnell (1953) ](https://learninghub.primer-e.com/link/954#bkmrk-horsnell1953)
Horsnell, G. (1953) The effect of unequal group variances on the *F*-test for the homogeneity of group means. *Biometrika*, **40**, 128-136.
---
[ Hotelling (1947) ](https://learninghub.primer-e.com/link/954#bkmrk-hotelling1947)
Hotelling, H. (1947) Multivariate quality control illustrated by air testing of sample bombsights. pp. 111-184 In: *Techniques of statistical analysis* (eds Eisenhart, C., Hastay, M.W. & Wallis, W.A.). New York, NY, USA: McGraw Hill.
---
[ Jensen *et al.* (2006) ](https://learninghub.primer-e.com/link/954#bkmrk-jensenetal2006)
Jensen, W.A., Jones-Farmer, L.A., Champ, C.W. & Woodall, W.H. (2006) Effects of parameter estimation on control chart properties: a literature review. *Journal of Quality Technology*, **38**, 349-364.
---
[ Johnson & Weerhandi (1988) ](https://learninghub.primer-e.com/link/954#bkmrk-johnsonweerhandi1988)
Johnson, R.A. & Weerhandi, S. (1988) A Bayesian solution to the multivariate Behrens-Fisher problem. *Journal of the American Statistical Association*, **83**, 145-149.
---
[ Kendall (1938) ](https://learninghub.primer-e.com/link/954#bkmrk-kendall1938)
Kendall, M.G. (1938) A new measure of rank correlation. *Biometrika*, **30**, 81-93.
---
[ Kendall (1945) ](https://learninghub.primer-e.com/link/954#bkmrk-kendall1945)
Kendall, M. G. (1945) The treatment of ties in ranking problems. *Biometrika*, **33**, 239-251.
---
[ Kolmogorov (1933) ](https://learninghub.primer-e.com/link/954#bkmrk-kolmogorov1933)
Kolmogorv, A. (1933) Sulla determinazione empirica di una legge di distribuzione. *Giornale dell'Instituto Italiano degli Attuari*, **4**, 83-91.
---
[ Kolmogorov (1941) ](https://learninghub.primer-e.com/link/954#bkmrk-kolmogorov1941)
Kolmogorv, A. (1941) Confidence limits for an unknown distribution function. *The Annals of Mathematical Statistics*, **12**, 461-463.
---
[ Krishnamoorthy & Lu (2010) ](https://learninghub.primer-e.com/link/954#bkmrk-krishnamoorthylu2010)
Krishnamoorthy, K. & Lu, F. (2010) A parametric bootstrap solution to the MANOVA under heteroscedasticity. *Journal of Statistical Computation and Simulation*, **80**, 873-997.
---
[ Kruskal (1952) ](https://learninghub.primer-e.com/link/954#bkmrk-kruskal1952)
Kruskal, W.H. (1952) A nonparametric test for the several sample problem. *The Annals of Mathematical Statistics*, **23**, 525-540.
---
[ Kruskal & Wallis (1952) ](https://learninghub.primer-e.com/link/954#bkmrk-kruskalwallis1952)
Kruskal, W.H. & Wallis, W.A. (1952) Use of ranks in one-criterion variance analysis. *Journal of the American Statistical Association*, **47**, 583-621.
---
[ Kruskal & Wish (1978) ](https://learninghub.primer-e.com/link/954#bkmrk-kruskalwish1978)
Kruskal, J.B. & Wish, M. (1978) *Multidimensional scaling*. Beverly Hills, CA, USA: Sage Publications.
---
[ Ledoit & Wolf (2003) ](https://learninghub.primer-e.com/link/954#bkmrk-ledoitwolf2003)
Ledoit, O. & Wolf, M. (2003) Improved estimation of the covariance matrix of stock returns with an application to portfolio selection. *Journal of Empirical Finance*, **10**, 603-621.
---
[ Ledoit & Wolf (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-ledoitwolf2004)
Ledoit, O. & Wolf, M. (2004) A well-conditioned estimator for large-dimensional covariance matrices. *Journal of Multivariate Analysis*, **88**, 365-411.
---
[ Mac Nally & Timewell (2005) ](https://learninghub.primer-e.com/link/954#bkmrk-macnallytimewell2005)
Mac Nally, R. & Timewell, C.A.R. (2005) Resource availability controls bird-assemblage composition through interspecific agression. *Auk*, **122**, 1097-1111.
---
[ Manly (2006) ](https://learninghub.primer-e.com/link/954#bkmrk-manly2006)
Manly, B.F.J. (2006) *Randomization, bootstrap and Monte Carlo methods in biology, 3rd edition*. London, UK: Chapman and Hall.
---
[ Mann & Whitney (1947) ](https://learninghub.primer-e.com/link/954#bkmrk-mannwhitney1947)
Mann, H.B. & Whitney, D.R. (1947) On a test of whether one of two random variables is stochastically larger than the other. *The Annals of Mathematical Statistics*, **18(1)**, 50-60.
---
[ McArdle & Anderson (2001) ](https://learninghub.primer-e.com/link/954#bkmrk-mcardleanderson2001)
McArdle B.H. & Anderson, M.J. (2001) Fitting multivariate models to community data: a comment on distance-based redundancy analysis. *Ecology*, **82**, 290-297.
---
[ McArdle & Anderson (2004) ](https://learninghub.primer-e.com/link/954#bkmrk-mcardleanderson2004)
McArdle B.H. & Anderson, M.J. (2004) Variance heterogeneity, transformations, and models of species abundance: a cautionary tale. *Canadian Journal of Fisheries and Aquatic Sciences*, **61**, 1294-1302.
---
[ Mead (1988) ](https://learninghub.primer-e.com/link/954#bkmrk-mead1988)
Mead, R. (1988) *The design of experiments: statistical principles for practical application*. Cambridge, UK: Cambridge University Press.
---
[ Montgomery (2020) ](https://learninghub.primer-e.com/link/954#bkmrk-montgomery2020)
Montgomery, D.C. (2020) *Introduction to statistical quality control, 8th edition*. New York, NY, USA: John Wiley & Sons.
---
[ Myers *et al.* (2021) ](https://learninghub.primer-e.com/link/954#bkmrk-myersetal2021)
Myers, E.M.V., Eme, D., Liggins, L., Harvey, E.S., Roberts, C.D. & Anderson, M.J. (2021) Functional beta diversity of New Zealand fishes: characterising morphological turnover along depth and latitude gradients, with derivation of functional bioregions. *Austral Ecology*, **46**, 965-981.
[link to source](https://onlinelibrary.wiley.com/doi/10.1111/aec.13078)
---
[ Opgen-Rhein & Strimmer (2007) ](https://learninghub.primer-e.com/link/954#bkmrk-opgenrheinstrimmer2007)
Opgen-Rhein, R. & Strimmer, K. (2007) Accurate ranking of differentially expressed genes by a distribution-free shrinkage approach. *Statistical Applications in Genetics and Molecular Biology*, **6(1)**, article 9.
---
[ Parzen (1962) ](https://learninghub.primer-e.com/link/954#bkmrk-parzen1962)
Parzen, E. (1962) On estimation of a probability density function and mode. *The Annals of Mathematical Statistics*, **33**, 1065-1076.
---
[ Pearson & Blackstock (1984) ](https://learninghub.primer-e.com/link/954#bkmrk-pearsonblackstock1984)
Pearson, T.H. & Blackstock, J. (1984) *Garroch Head sludge dumping ground survey, final report*. Dunstaffnage Marine Research Laboratory, UK (unpublished).
---
[ Pratt (1959) ](https://learninghub.primer-e.com/link/954#bkmrk-pratt1959)
Pratt, J.W. (1959) Remarks on zeros and ties in the Wilcoxon signed rank procedures. *Journal of the American Statistical Association*, **54**, 655-667.
---
[ Prober *et al.* (2007) ](https://learninghub.primer-e.com/link/954#bkmrk-proberetal2007)
Prober, S.M., Thiele, K.R. & Lunt, I.D. (2007) Fire frequency regulates tussock grass composition, structure and resilience in endangered temperate woodlands. *Austral Ecology*, **32**, 808-824.
---
[ Quesenberry (2007) ](https://learninghub.primer-e.com/link/954#bkmrk-quesenberry2007)
Quesenberry, C.P. (2007) *SPC methods for quality improvement*. New York, NY, USA: John Wiley & Sons.
---
[ Rao (1968) ](https://learninghub.primer-e.com/link/954#bkmrk-rao1968)
Rao, J.N.K. (1968) On expectations, variances, and covariances of ANOVA mean squares by 'synthesis'. *Biometrics*, **24**, 963-978.
---
[ Rosenblatt (1956) ](https://learninghub.primer-e.com/link/954#bkmrk-rosenblatt1956)
Rosenblatt, M. (1956) Remarks on some nonparametric estimates of a density function. *The Annals of Mathematical Statistics*, **27**, 832-837.
---
[ Sammon (1969) ](https://learninghub.primer-e.com/link/954#bkmrk-sammon1969)
Sammon, J.W. (1969) A nonlinear mapping for data structure analysis. *IEEE Transactions on Computers*, **18**, 401-409.
---
[ Satterthwaite (1941) ](https://learninghub.primer-e.com/link/954#bkmrk-satterthwaite1941)
Satterthwaite, F.E. (1941) Synthesis of variance. *Psychometrika*, **6**, 309-316.
---
[ Sauer *et al.* (2019) ](https://learninghub.primer-e.com/link/954#bkmrk-saueretal2019)
Sauer, J.R., Niven, D.K., Hines, J.E., Ziolkowski, D.J., Jr., Pardieck, K.L., Fallon, J.E. & Link, W.A. (2019) *The North American Breeding Bird Survey, Results and Analysis 1966 - 2019 (Version 2.07.2019)*. USGS Patuxent Wildlife Research Center, Laurel, MD, USA.
---
[ Schäfer & Strimmer (2005) ](https://learninghub.primer-e.com/link/954#bkmrk-schaferstrimmer2005)
Schäfer, J. & Strimmer, K. (2005) A shrinkage approach to large-scale covariance matrix estimation and implications for functional genomics. *Statistical Applications in Genetics and Molecular Biology*, **4(1)**, article 32.
---
[ Schenk *et al.* (2025) ](https://learninghub.primer-e.com/link/954#bkmrk-schenketal2025)
Schenk, S., Supratya, V.P., Martone, P.T. & Parfrey, L.W. (2025) Monthly macroalgal surveys reveal a diverse and dynamic community in an urban intertidal zone. *Botany*, **103**, 1-12.
---
[ Seber (1984) ](https://learninghub.primer-e.com/link/954#bkmrk-seber1984)
Seber, G.A.F. (1984) *Multivariate observations*. New York, NY, USA: John Wiley & Sons.
---
[ Sheskin (2011) ](https://learninghub.primer-e.com/link/954#bkmrk-sheskin2011)
Sheskin, D.J. (2011) *Handbook of parametric and nonparametric statistical procedures, 5th edition*. Boca Raton, FL, USA: Chapman & Hall/CRC Press, Taylor & Francis Group.
---
[ Shewhart (1931) ](https://learninghub.primer-e.com/link/954#bkmrk-shewhart1931)
Shewhart, W.A. (1931) *Economic control of quality of manufactured product*. New York, NY, USA: Van Nostrand.
---
[ Shewhart (1939) ](https://learninghub.primer-e.com/link/954#bkmrk-shewhart1939)
Shewhart, W.A. (1939) *Statistical method from the viewpoint of quality control*. Graduate School of the Department of Agriculture, Washington, D.C., USA.
---
[ Silverman (1986) ](https://learninghub.primer-e.com/link/954#bkmrk-silverman1986)
Silverman, B.W. (1986) *Density estimation for statistics and data analysis*. London, UK: Chapman & Hall.
---
[ Smirnov (1939a) ](https://learninghub.primer-e.com/link/954#bkmrk-smirnov1939a)
Smirnov, N.V. (1939a) On the estimation of the discrepancy between empirical curves of distribution for two independent samples. *Bulletin de l'Université de Moscou, Série internationale (Mathématiques)*, **2(2)**, 16 pp.
---
[ Smirnov (1939b) ](https://learninghub.primer-e.com/link/954#bkmrk-smirnov1939b)
Smirnov, N.V. (1939b) Sur les écarts de la courbe de distribution empirique. *Recueil Mathématiques de Moscou (Matematicheskii Sbornik)*, **6(48)**, 3-26.
---
[ Smith & Anderson (2016) ](https://learninghub.primer-e.com/link/954#bkmrk-smithanderson2016)
Smith, A.N.H. & Anderson, M.J. (2016) Marine reserves indirectly affect fine-scale habitat associations, but not overall densities, of small benthic fishes. *Ecology and Evolution*, **6**, 6648-6661.
---
[ Smith *et al.* (2014) ](https://learninghub.primer-e.com/link/954#bkmrk-smithetal2014)
Smith, A.N.H., Anderson, M.J., Millar, R.B. & Willis, T.J. (2014) Effects of marine reserves in the context of spatial and temporal variation: an analysis using Bayesian zero-inflated mixed models. *Marine Ecology Progress Series*, **499**, 203-216.
---
[ Snedecor (1946) ](https://learninghub.primer-e.com/link/954#bkmrk-snedecor1946)
Snedecor, G.W. (1946) *Statistical methods, 4th edition*. Ames, Iowa, USA: Iowa State College Press.
---
[ Somerfield *et al.* (1994a) ](https://learninghub.primer-e.com/link/954#bkmrk-somerfieldetal1994a)
Somerfield, P.J., Gee, J.M. & Warwick, R.M. (1994a) Benthic community structure in relation to an instantaneous discharge of waste water from a tin mine. *Marine Pollution Bulletin*, **28**, 363-369.
---
[ Somerfield *et al.* (1994b) ](https://learninghub.primer-e.com/link/954#bkmrk-somerfieldetal1994b)
Somerfield, P.J., Gee, J.M. & Warwick, R.M. (1994b) Soft sediment meiofaunal community structure in relation to a long-term heavy metal gradient in the Fal estuary system. *Marine Ecology Progress Series*, **105**, 79-88.
---
[ Somerfield & Clarke (2013) ](https://learninghub.primer-e.com/link/954#bkmrk-somerfieldclarke2013)
Somerfield, P.J. & Clarke, K.R. (2013) Inverse analysis in non-parametric multivariate analyses: distinguishing groups of associated species which covary coherently across samples. *Journal of Experimental Marine Ecology and Biology*, **449**, 261-273.
---
[ Somerfield *et al.* (2021a) ](https://learninghub.primer-e.com/link/954#bkmrk-somerfieldetal2021a)
Somerfield, P.J., Clarke, K.R. & Gorley, R.N. (2021a) A generalised analysis of similarities (ANOSM) statistic for designs with ordered factors. *Austral Ecology*, **46**, 901-910.
---
[ Somerfield *et al.* (2021b) ](https://learninghub.primer-e.com/link/954#bkmrk-somerfieldetal2021b)
Somerfield, P.J., Clarke, K.R. & Gorley, R.N. (2021a) Analysis of similarities (ANOSM) for 2-way layouts using a generalised ANOSIM statistic, with comparative notes on Permutational Multivariate Analysis of Variance (PERMANOVA). *Austral Ecology*, **46**, 911-926.
---
[ Trott (2022) ](https://learninghub.primer-e.com/link/954#bkmrk-trott2022)
Trott, T.J. (2022) Mesoscale spatial patterns of Gulf of Maine rocky intertidal communities. *Diversity*, **14**, 557.
[link to source](https://doi.org/10.3390/d14070557)
---
[ Terlizzi *et al.* (2005) ](https://learninghub.primer-e.com/link/954#bkmrk-terlizzietal2005)
Terlizzi, A., Scuderi, D., Fraschetti, S. & Anderson, M.J. (2005) Quantifying effects of pollution on biodiversity: a case study of highly diverse molluscan assemblages in the Mediterranean. *Marine Biology*, **148**, 293-305.
---
[ Ullah *et al.* (2017) ](https://learninghub.primer-e.com/link/954#bkmrk-ullahetal2017)
Ullah, I., Pawley, M.D.P., Smith, A.N.H. & Jones, B. (2017) Improving the detection of unusual observations in high-dimensional settings. *Australian & New Zealand Journal of Statistics*, **59**, 449-462.
---
[ Underwood (1991) ](https://learninghub.primer-e.com/link/954#bkmrk-underwood1991)
Underwood, A.J. (1991) Beyond BACI: experimental designs for detecting human environmental impacts on temporal variations in natural populations. *Australian Journal of Marine and Freshwater Research*, **42**, 569-587.
---
[ Underwood (1992) ](https://learninghub.primer-e.com/link/954#bkmrk-underwood1992)
Underwood, A.J. (1992) Beyond BACI: the detection of environmental impcts on populations in the real, but variable, world. *Journal of Experimental Marine Biology and Ecology*, **161**, 145-178.
---
[ Underwood (1994) ](https://learninghub.primer-e.com/link/954#bkmrk-underwood1994)
Underwood, A.J. (1994) On beyond BACI: sampling designs that might reliably detect environmental disturbances. *Ecological Applications*, **4**, 3-15.
---
[ Underwood (1997) ](https://learninghub.primer-e.com/link/954#bkmrk-underwood1997)
Underwood, A.J. (1997) *Experiments in ecology: their logical design and interpretation using analysis of variance*. Cambridge, UK: Cambridge University Press.
---
[ Veale *et al.* (2014) ](https://learninghub.primer-e.com/link/954#bkmrk-vealeetal2014)
Veale, L., Tweedley, J.R., Clarke, K.R., Hallett, C.S. & Potter, I.C. (2014) Characteristics of the ichthyofauna of a temperate microtidal estuary with a reverse salinity gradient, including inter-decadal comparisons. *Journal of Fish Biology*, **85**, 1320-1354.
---
[ Wang (1971) ](https://learninghub.primer-e.com/link/954#bkmrk-wang1971)
Wang, Y.Y. (1971) Probabilities of the type I errors of the Welch tests for the Behrens–Fisher problem. *Journal of the American Statistical Association*, **66**, 605-608.
---
[ Warwick & Clarke (1993) ](https://learninghub.primer-e.com/link/954#bkmrk-warwickclarke1993)
Warwick, R.M. & Clarke, K.R. (1993) Increased variability as a symptom of stress in marine communities. *Journal of Experimental Marine Biology and Ecology*, **172**, 215-226.
---
[ Weerhandi (1993) ](https://learninghub.primer-e.com/link/954#bkmrk-weerhandi1993)
Weerhandi, S. (1993) Generalized confidence intervals. *Journal of the American Statistical Association*, **88**, 899-905.
---
[ Welch (1938) ](https://learninghub.primer-e.com/link/954#bkmrk-welch1938)
Welch, B.L. (1938) The significance of the difference between two means when population variances are unequal. *Biometrika*, **29**, 350-362.
---
[ Whittaker (1952) ](https://learninghub.primer-e.com/link/954#bkmrk-whittaker1952)
Whittaker, R.H. (1952) A study of summer foliage insect communities in the Great Smoky Mountains. *Ecological Monographs*, **22**, 1-44.
---
[ Wilcoxon (1945) ](https://learninghub.primer-e.com/link/954#bkmrk-wilcoxon1945)
Wilcoxon, F. (1945) Individual comparisons by ranking methods. *Biometrics Bulletin*, **1(6)**, 80-83.
---
[ Wilcoxon (1949) ](https://learninghub.primer-e.com/link/954#bkmrk-wilcoxon1949)
Wilcoxon, F. (1949) *Some rapid approximate statistical procedures*. New York, NY, USA: American Cyanamid Corporation.
---
[ Winsor & Clarke (1940) ](https://learninghub.primer-e.com/link/954#bkmrk-winsorclarke1940)
Winsor, C.P. & Clarke, G.L. (1940) A statistical study of variation in the catch of plankton nets. *Journal of Marine Research*, **3**, 1–34.
---
[ Wong (2011) ](https://learninghub.primer-e.com/link/954#bkmrk-wong2011)
Wong, B. (2011) Points of view: colorblindness. *Nature Methods*, **8**, 441.
---