Background: Phenotypic Similarity

To soften the Jaccard and Chao-Jaccard indices such that they include additional information about the resemblance of two samples, we first need a separate measure that can be used to establish the morphological similarity of a given species in two different samples. For example, if we know that an arbitrary species x exists in sample A and B but its morphology varies (for example, x could be twice as large in B), incorporating that information using the soft similarity approach described necessitates an appropriate measure of morphological similarity. Assume individual taxa can be represented as sets C and D consisting of k morphological measurements (C = {c ₁, c ₂, …c _k} and D = {d ₁, d ₂, …d _k}), where c _k and d _k are the mean values of trait k for that species in that sample. Here we are primarily concerned with metric morphological measurements (e.g., body size/volume) and thus primarily use the term “morphology” when describing the set of measurements used for each taxon, but the k measurements could be any sort of phenotypic measurement (fractional trophic level, functional ecospace characters [sensu Novack-Gottshall Reference Novack-Gottshall2007], etc.). Because we are allowing for mixed data types when describing the morphology/phenotype of a given species, a measure that can handle different types of data is needed. Gower (Reference Gower1971) generalized Jaccard's concept of set-based similarity to include the ability to compare collections of objects/traits in which the types of elements in a set could vary. Given that phenotypic measurements may vary in scale and in the type of measurement (e.g., metric, binary, ordinal), we utilize the Gower index as an approach to quantifying phenotypic similarity. Gower (Reference Gower1971; see also Podani Reference Podani1999) proposed a suite of sub-indices for comparing different types of characters/elements of a given set (e.g., binary, nominal, numeric, and ordinal values), which would be summed to produce the final similarity index. Given this, we can calculate similarity as:

(1)$$ {Gower( {C, \;D} ) = \displaystyle{{\vert {C\,\tilde{\cap }\,D} \vert } \over {\vert {C_\delta \cap D_\delta } \vert }}.\;} $$

Modifying the general argument described in Gower (Reference Gower1971), the similarity of species C and D is calculated using a membership function that assigns each element c _k a similarity value to its counterpart d _k, based on the class of element being compared:

(2)$$ {\vert {C\,\tilde{\cap }\,D} \vert = \mathop \sum \limits_{k = 1}^{\vert C \vert } \mu ( {c_k, \;C\,\tilde{\cap }D}\, ) \;{\rm for}\;{\rm all}\;c_k\in C.} $$

In the original formulation of the Gower index, the membership of a given character c _k in the intersection of C and D was calculated using three base cases:

(3)$${{\rm \mu }( {c_k, \;C\,\tilde{\cap }\,D}) = \left\{{\matrix{ {{\rm \mu }_{Interval}( {c_k, \;C\,\tilde{\cap }\,D} ) \,{\rm if}\;c_k\;{\rm is}\;\,{\rm on}\;\,{\rm the}\,{\rm \;interval}\;\,{\rm or}\,{\rm \;ratio}\,\,{\rm scale}} \cr \hskip-8.5pc{{\rm \mu }_{Binary}( {c_k, \;C\,\tilde{\cap }\,D} ) \;\,{\rm if}\;c_k\;{\rm is}\,{\rm \;binary}} \cr\hskip-7.5pc {{\rm \mu }_{Nominal}( {c_k, \;C\,\tilde{\cap }\,D} ) \;\,{\rm if}\;c_k\;{\rm is}\,\;{\rm nominal}} \cr } } \right.\;.}$$

In the case of numeric data (e.g., intervals, ratios, ordinal characters, real numbers), the similarity between two traits c _k and d _k was calculated as:

(4)$$ {{\rm \mu }_{Interval}( {c_k, \;C\,\tilde{\cap }\,D} ) = 1-\displaystyle{{\vert {\vert {c_k-d_k} \vert } \vert } \over {R( k ) }}, \;\;} $$

where R(k) is the range of values observed for trait k of the set. Note that there are important features of the Gower index that must be accounted for if it is used as a measurement of phenotypic similarity. The first is the problem of interspecific variation. Most existing implementations of the Gower index use the range of all values for a given trait to calculate the similarity of any two individual trait values. For example, if all species being compared were measured in terms of maximum linear dimension (as a proxy for body size), then the proximity of any two body-size measurements for the same species is taken in relation to the range of body sizes for all species. What this means is that the relative proximity of any two measurements is driven by the absolute difference between the smallest and largest species in the two samples. Consider, for example, the proximity of two arbitrarily chosen mean body-size measurements (c _k = 1.2 mm and d _k = 2.1 mm). If, as stated above, the value of trait k is the mean maximum linear dimension for that species, consider the following two cases. In case 1, the smallest maximum linear dimension observed in the pooled assemblage is min _k = 1 mm, while the largest maximum linear dimension observed in the pooled assemblage is max _k = 5 mm. Here we scale the values between zero and one by dividing by the maximum value. In this case, the proximity of c _k and d _k is 0.78 ($1-{{\vert {\vert {0.24-0.42} \vert } \vert } \over {0.8}}$). Now consider case 2, where the largest maximum linear dimension value in the pooled assemblage is 10 times greater (max _k = 50 mm). In this case, the proximity of c _k and d _k is 0.98 ($1-{{\vert {\vert {0.024-0.042} \vert } \vert } \over {0.98}}$). This general pattern shows how the relative proximity (and contribution to the final similarity value) of any two measurements increases with the total range of values and decreases as the range decreases. This makes interpreting the proximity of any two values difficult—is the range of interspecies body sizes in each sample a meaningful space to compare the proximity of intraspecies body sizes? We suggest that when quantifying the similarity of a given species to itself in another sample, it is not, and is particularly problematic when interspecies measurements differ by orders of magnitude. For example, if a species of zooplankton is twice as large in one sample, but the proximity of these two measurements is compared on a scale that incorporates body-size measurements for something several orders of magnitude larger (e.g., a fish), the phenotype of the zooplankton will be treated as if it was essentially exactly the same (even though it doubled/halved, depending on your frame of reference). As such, we suggest using the range of values observed for that species to compare the proximity of mean values for that species. In these cases, the absolute difference between the mean values for a given trait for each species in each sample is divided by the absolute difference between the minimum and maximum values observed for that species in the pooled assemblage. Note that this suggestion is made in relation to calculating the new indices described here only—in other cases (e.g., measuring functional diversity; Leinster and Cobbold Reference Leinster and Cobbold2012), the proximity of measurements should be made on the interspecific scale.

In addition to this problem of scaling, there exists an issue related to handling binary/logical/nominal characters when they vary in their expression between individuals within a given sample (e.g., a feature like a marking could be present on some individuals but not others). While the notion of a mean value for a metric phenotypic character (that is, one that is expressed as a real number/integer) is readily understood as the usual meaning of the term, it is less clear how we should handle binary/logical/nominal traits. In the standard Gower index, conventionally used to measure the dissimilarity between two phenotypes where the trait has a fixed value of either true/false or present/absent, the case for binary variables is:

(5)$$ {{\rm \mu }_{Binary}( {c_k, \;C\,\tilde{\cap }\,D} ) = \left\{{\matrix{ {1\;\,{\rm if}\;c_k = d_k = 1} \cr {0\;\,{\rm if}\;c_k = 1\;\,{\rm and}\,\;d_k = 0} \cr {0\;\,{\rm if}\,c_k = 0\;\,{\rm and}\,\;d_k = 1} \cr {0\,{\rm \;if}\,\;c_k = d_k = 0} \cr } } \right..} $$

Finally, the case for nominal variables is:

(6)$${{\rm \mu }_{Nominal}( {c_k, \;C\,\tilde{\cap }\,D} ) = \left\{{\matrix{ {1\,{\rm \;if}\,\;c_k = d_k} \cr {0\;\,{\rm if}\;\,c_k\ne d_{k.}} \cr } } \right..} $$

When used to describe a given species, Gower's approach to binary/nominal characters assumes they are a fixed component of that species, found on all individuals. When comparing the mean phenotype of a species in two different samples, we suggest converting the binary/logical/nominal characters into frequencies (the probability that a randomly chosen individual of a given species in one sample has the feature) and comparing the frequencies using the case for interval data described above. Start by recording the presence/absence of that binary trait for each individual. Next, divide the total number of positive cases (“present” values) by the total number of individuals of that species, to attain a frequency for the binary trait. Next, calculate the similarity of the frequency of that trait in each sample using the $\mu _{Interval}( {c_k, \;C\tilde{\cap }D} ) $ formula described above. In these cases, if the frequency of the binary trait is equal in both samples, it contributes a value of 1. If the frequency of the binary trait is 1 in one sample and 0 in another, the trait contributes a value of 0. In all other cases, it contributes a value between 0 and 1 equal to the similarity of the frequency that the trait occurs at in the two samples. In a practical sense, binary/logical/unordered nominal characters should be converted into columns of 1s and 0s, where 1 indicates the character is present on that individual (i.e., it is a 1 for a binary character, TRUE for a logical character, and present for a nominal character). If we consider only the values assigned to individuals of a single species, the mean value of these columns will be the frequency of true positives (the probability that a randomly selected individual of that species has the trait). The included R code (see Supplementary Material and the Appendix) automatically calculates this modified version of the Gower similarity based on the type of data and can handle binary, logical, and nominal traits that vary between individual specimens.

Soft Jaccard Similarity

The classic Jaccard index uses a strict cardinality when counting the number of elements in the intersection—the cardinality of the intersection increases by 1 each time a shared species is identified (Jimenez et al. Reference Jimenez, Gonzalez and Gelbukh2016) (Fig. 1). The soft Jaccard index proposed here uses a soft cardinality, wherein shared species increase the cardinality of the intersection by a value equal to their Gower similarity—if the phenotype of an arbitrary shared species pair x is identical in the two samples, they contribute a value of 1 (as in the strict cardinality) (Fig. 1). If their phenotype shows no overlap at all, they contribute a value of 0. When their phenotype is similar in the two samples, they contribute a value equal to their Gower similarity, which is expressed on the interval [0,1] (Fig. 1). Note that we are comparing the phenotype of a given species in one sample to the phenotype of the same species in another sample—we are not computing pairwise dissimilarities between all taxa in a sample as is frequently done in measurements of trait/functional diversity (see Leinster and Cobbold [Reference Leinster and Cobbold2012] for an example of the latter type of calculation). Repeating this process for every species that is found in both samples, we can derive a Jaccard-type index that weights incidence by morphological/phenotypic similarity. Imagine two or more samples taken from communities where n species can be described in terms of metrics associated with morphology/phenotype (height and thickness of sponges/corals, leaf area of plants, body volume, diet proportions, etc.). Begin with two sets of m shared species with associated phenotypic measurements F = {C ₁, C ₂, …C _m} and G = {D ₁, D ₂, …D _m}. Assume each shared species (denoted C _m in F and D _m in G) has a set of k phenotypic measurements associated with it (that is, C _m = {c ₁, c ₂, …c _k} and D _m = {d ₁, d ₂, …d _k}), where each measurement c _k and d _k can be expressed as above. Note that the characters used to measure one species do not need to be the same as the characters used to measure another. Our proposed index takes the form of the original Jaccard index, with the numerator replaced by a value equal to the soft cardinality of the intersection of sets F and G:

(7)$${Soft\;Jaccard( {F, \;G} ) = \displaystyle{{\vert {F\,\tilde{\cap }\,G} \vert } \over {\vert {F\,\tilde{\cup }\,G} \vert }}.\;} $$

The soft cardinality of the intersection of F and G is calculated using a membership function that weights the incidence of a given species m by its phenotypic dissimilarity between the two sites (that is, it softens the definition of “shared”):

(8)$$ {\vert {F\,\tilde{\cap }\,G} \vert = \mathop \sum \limits_{m = 1}^{\vert F \vert } {\rm \mu }( {C_m, \;F\,\tilde{\cap }\,G} ) \;\,{\rm for}\,{\rm \;all}\,\;C_m\in F.\;} $$

The membership of each shared species C _m in F in the soft intersection of F and G ($F\,\widetilde{\cap }\,G$) is a value equal to the phenotypic similarity of that species when compared with its counterpart in G (D _m):

(9)$$\mu \left( {C_m,F\,\widetilde{\cap}\,G} \right) = \displaystyle{{\left| {C_m\,\widetilde{\cap}\,D_m} \right|} \over {\left| {C_{m_\delta} \cap D_{m_\delta} } \right|}},$$

where phenotypic similarity is the Gower similarity of the two taxa, calculated as above. We divide this value by the cardinality of the soft union of F and G ($F\widetilde{\cup }G$), calculated as the number of species in F and G minus the species shared by the two samples, where the definition of shared is softened as above:

(10)$${\vert {F\,\tilde{\cup }\,G} \vert = \vert F \vert + \vert G \vert -\vert {F\,\tilde{\cap }\,G} \vert.\;} $$

Note that if the phenotype/morphology of all shared species is the same in both samples, the index reduces to the strict Jaccard index.

Figure 1.

Image showing how the cardinality of the intersection is calculated in the strict Jaccard index and the soft Jaccard index. Arrows and their associated values indicate the contribution of that species to the cardinality of the intersection. In the strict Jaccard index (left), shared species contribute a value of 1, regardless of any degree of morphological/phenotypic variation that might exist. In the soft Jaccard index (right), shared species contribute a value on the interval [0,1], depending on their morphological similarity. In this way, the soft Jaccard index incorporates more information than its strict counterpart. $\vert {A\cap B} \vert $ indicates the strict cardinality of the intersection, while $\vert {A\widetilde{\cap }B} \vert $ indicates the soft cardinality of the intersection. Values chosen for the soft Jaccard example are arbitrary and used for illustration purposes only.

Example

Assume X = {x ₁, x ₂, x ₃} and Y = {y ₁, y ₂, y ₃} are communities composed of the same three species (that is, x ₁ = y ₁, x ₂ = y ₂, and x ₃ = y ₃, such that X = Y). By definition $\vert {X\cap Y} \vert = 3$. We can represent the intersection of X and Y as an identity matrix I (wherein the diagonal indicates it is TRUE that species 1 is present in X and Y, etc.):

(11)$$I = \left[{\matrix{ 1 & 0 & 0 \cr 0 & 1 & 0 \cr 0 & 0 & 1 \cr } } \right], \;$$

such that $\sum ( I ) = \vert {X\cap Y} \vert $. When calculating the strict Jaccard similarity of X and Y, we take the cardinality of the intersection and divide it by the cardinality of the union ($\vert {X\cup Y} \vert = \vert X \vert + \vert Y \vert -\vert {X\cap Y} \vert $), such that $Jaccard( {X, \;Y} ) = {{\vert {X\cap Y} \vert } \over {\vert {X\cup Y} \vert }} = {3 \over 3} = 1$. To calculate the soft Jaccard similarity of X and Y, we replace the values in the diagonal of the identity matrix I with values equal to the phenotypic similarity of that species as quantified between X and Y (calculated as Gower(x _i, y _i)), producing a special diagonal matrix we are calling a soft identity matrix $\tilde{I}$:

(12)$${\tilde{I} = \left[{\matrix{ {Gower( {x_1, \;y_1} ) } & 0 & 0 \cr 0 & {Gower( {x_2, \;y_2} ) } & 0 \cr 0 & 0 & {Gower( {x_3, \;y_3} ) } \cr } } \right].}$$

For example, if Gower(x ₁, y ₁) = 0.97, Gower(x ₂, y ₂) = 0.33, and Gower(x ₃, y ₃) = 0.72, the soft identity matrix is:

(13)$${\tilde{I} = \left[{\matrix{ {0.97} & 0 & 0 \cr 0 & {0.33} & 0 \cr 0 & 0 & {0.72} \cr } } \right].}$$

We can calculate the soft cardinality of the intersection of X and Y as $\sum ( {\tilde{I}} ) = \vert {X\,\widetilde{\cap }\,Y} \vert $. Following from this, the soft cardinality of the union of X and Y is $\vert {X\,\widetilde{\cup }\,Y} \vert = \vert X \vert + \vert Y \vert -\vert {X\,\widetilde{\cap }\,Y} \vert $. Finally, the soft Jaccard similarity of X and Y is $Soft\;Jaccard( {X, \;Y} ) = {{\vert {X\,\tilde{\cap }\,Y} \vert } \over {\vert {X\,\tilde{\cup }\,Y} \vert }} = {{2.02} \over {3.98}} = 0.51$. Note that the use of the Gower index to calculate the similarity of shared elements is only one way to produce a soft identity matrix (chosen here for its applicability to questions of phenotypic similarity). A generalized form is:

(14)$${\tilde{I} = \left[{\matrix{ {S( {x_1, \;y_1} ) } & 0 & 0 \cr 0 & {S( {x_2, \;y_2} ) } & 0 \cr 0 & 0 & {S( {x_3, \;y_3} ) } \cr } } \right], \;} $$

where S(x _i, y _i) can be any similarity index expressed on the interval [0,1] (e.g., see Molodtsov Reference Molodtsov1999; Chum et al. Reference Chum, Philbin and Zisserman2008, Jimenez et al. Reference Jimenez, Gonzalez and Gelbukh2016; Chauvin et al. Reference Chauvin, Kumar, Wachinger, Vangel, de Guise, Desrosiers, Wells and Toews2020).

Soft Chao-Jaccard Similarity

Chao et al. (Reference Chao, Chazdon, Colwell and Shen2005, Reference Chao, Chazdon, Colwell and Shen2006) noted that classic incidence-based indices like Jaccard and Sørenson-Dice are sensitive to sample size and the presence of rare species. Because of this, rare species and common species are treated equally, despite sound a priori reasons to compare samples in terms of individuals, rather than species/clades (Magurran Reference Magurran2004). This lack of appropriate feature weighting can lead the classic Jaccard index to be downward biased when comparing samples taken from populations with rare species (Wolda Reference Wolda1981). Chao et al.'s (Reference Chao, Chazdon, Colwell and Shen2005) index circumvents this issue by working with sets U = {u ₁, u ₂, …u _n} and V = {v ₁, v ₂, …v _n} of relative abundances u _n and v _n of n shared species. Here we describe a method for producing a soft Chao-Jaccard index, wherein variation in abundance and intraspecific morphology is taken into account when comparing samples. As above, we imagine two or more samples taken from communities where n species can be described in terms of metrics associated with phenotype. If each species in each sample also has an associated abundance value, we can derive a soft version of the Chao-Jaccard similarity index. Given two sets F and G of m shared species with associated phenotypic measurements F = {C ₁, C ₂, …C _m} and G = {D ₁, D ₂, …D _m} and sample relative abundances U = {u ₁, u ₂, …u _m} and V = {v ₁, v ₂, …v _m}, with a set W of corresponding phenotypic similarity values for each shared species:

(15)$$\eqalign{W &= { {Gower( {C_1, \;D_1} ) , \;Gower( {C_2, \;D_2} ) , }}\cr & {\ldots Gower( {C_m, \;D_m} ) } , }\;$$

We can derive a soft version of the Chao-Jaccard index. Note that U and V are sets of relative abundances u _n and v _n of n shared species. Specifically, if P and Q are sets of absolute abundance values (specimen counts) associated with the n species found in U and V, respectively (P = {p ₁, p ₂, …, p _n}, Q = {q ₁, q ₂, …, q _n}), and a _n and b _n are presence/absence values associated with species n (such that a _n indicates whether species n is present in P [p _n > 0] and b _n indicates whether that species is found in Q [q _n > 0]), then $u_n = {{p_n} \over {\sum ( P ) }}( {a_n \times b_n} ) $. That is, u _n is equal to the relative abundance of species n in sample P if species n is a shared species, and equal to 0 if species n is not a shared species. The values included in set V are calculated in the same way with appropriate substitutions. To soften the original Chao-Jaccard index, we first weight the probability of selecting an individual representing a shared species from each sample using the phenotypic similarity of that species to its counterpart in the other sample. In doing this we are producing a soft cardinality of the sets U and V, denoted $\tilde{U}$ and $\tilde{V}$ and calculated as $\tilde{U} = U \times W$ and $\tilde{V} = V \times W$ respectively. From this, we insert the weighted relative abundances into the original Chao-Jaccard index to produce a measure of similarity sensitive to differences in abundance and mean phenotype of shared species:

(16)$$Soft\;Chao-Jaccard( {U, \;V} ) = \displaystyle{{\sum ( {\tilde{U}} ) \times \sum ( {\tilde{V}} ) } \over {\sum ( {\tilde{U}} ) + \sum ( {\tilde{V}} ) -\left({\sum ( {\tilde{U}} ) \times \sum ( {\tilde{V}} ) } \right)}}.$$

While the soft Jaccard index uses a soft cardinality of the intersection, the soft Chao-Jaccard index uses a soft cardinality of each individual sample. In essence, when counting the number of individuals of a given species in a given sample, we take into account the extent to which that species has a similar morphology in both samples. As above, this comparison is only taking into account the phenotype of the same species in the other samples. If all species have a maximum pairwise similarity of 1, the index reduces to the standard Chao-Jaccard index.

Example

As above, assume X = {x ₁, x ₂, x ₃} and Y = {y ₁, y ₂, y ₃} are communities composed of the same three species (that is, x ₁ = y ₁, x ₂ = y ₂, and x ₃ = y ₃ such that X = Y). Here we assign each species a relative abundance for that sample, and, as above, use U and V to indicate sets of relative abundances u _i and v _i. As we did with the intersection, we can represent U and V as diagonal matrices J and K (here the diagonal indicates the relative abundance of species 1 in community X, etc.). We use J and K to indicate diagonal matrices referring to individual sets, to differentiate them from the symbol I used above to represent the identity matrix for the intersection. For example:

(17)$$J = \left[{\matrix{ {u_1} & 0 & 0 \cr 0 & {u_2} & 0 \cr 0 & 0 & {u_3} \cr } } \right]\;\,{\rm and}\;K = \left[{\matrix{ {v_1} & 0 & 0 \cr 0 & {v_2} & 0 \cr 0 & 0 & {v_3} \cr } } \right]. $$

Note that as stated in Chao et al. (Reference Chao, Chazdon, Colwell and Shen2005) and described above, U and V include the relative abundances of shared species only—thus, there are no zeros on the diagonal of J and K. For sake of demonstration, assume u ₁ = v ₁ = 0.2, u ₂ = v ₂ = 0.2, and u ₃ = v ₃ = 0.6:

(18)$$J = K = \left[{\matrix{ {0.2} & 0 & 0 \cr 0 & {0.2} & 0 \cr 0 & 0 & {0.6} \cr } } \right]. $$

As was done above when calculating the soft Jaccard similarity, we first calculate a soft identity matrix for the intersection of X and Y (restating W as a diagonal matrix):

(19)$$\tilde{I} = \left[{\matrix{ {Gower( {x_1, \;y_1} ) } & 0 & 0 \cr 0 & {Gower( {x_2, \;y_2} ) } & 0 \cr 0 & 0 & {Gower( {x_3, \;y_3} ) } \cr } } \right].$$

For sake of demonstration, assume, as above, that Gower(x ₁, y ₁) = 0.97, Gower(x ₂, y ₂) = 0.33, and Gower(x ₃, y ₃) = 0.72, such that:

(20)$$\tilde{I} = \left[{\matrix{ {0.97} & 0 & 0 \cr 0 & {0.33} & 0 \cr 0 & 0 & {0.72} \cr } } \right].$$

We calculate the soft cardinality of U ($\vert {\tilde{U}} \vert $) and V ($\vert {\tilde{V}} \vert $) as follows:

(21)$$\tilde{U} = \tilde{I} \times J = \left[{\matrix{ {0.97} & 0 & 0 \cr 0 & {0.33} & 0 \cr 0 & 0 & {0.72} \cr } } \right]\times \left[{\matrix{ {0.2} & 0 & 0 \cr 0 & {0.2} & 0 \cr 0 & 0 & {0.6} \cr } } \right] = \left[{\matrix{ {0.194} & 0 & 0 \cr 0 & {0.066} & 0 \cr 0 & 0 & {0.432} \cr } } \right], $$

given that J = K, $\tilde{U} = \tilde{V}$. As above, $\vert {\tilde{U}} \vert = \sum ( {\tilde{U}} ) $ and $\vert {\tilde{V}} \vert = \sum ( {\tilde{V}} ) $. Following from this, the soft Chao-Jaccard similarity of X and Y is Soft Chao-Jaccard (U, V) = $ {{\sum ( {\tilde{U}} ) \times \sum ( {\tilde{V}} ) } \over {\sum ( {\tilde{U}} ) + \sum ( {\tilde{V}} ) -\left({\sum ( {\tilde{U}} ) \times \sum ( {\tilde{V}} ) } \right)}}=$ $ {{0.479} \over {0.905}} = 0.53$.

The soft Chao-Jaccard similarity of these two sets is slightly higher than the soft Jaccard similarity, because the most dissimilar species (species 2) is less abundant, and therefore contributes less to the final value.

Description of Datasets

To explore the new approach, we compared principal coordinates analysis (PCoA) outputs based on the original Chao-Jaccard index (Chao et al. Reference Chao, Chazdon, Colwell and Shen2005) with ones based on the soft Chao-Jaccard index (proposed here). Note that although we are choosing the Chao-Jaccard index for comparative purposes, we must emphasize that increasing the information provided to any index improves a researcher's ability to recognize short gradients using ordination techniques and does not reflect issues with the original index. Random samples were compared from two settings—a simulated river velocity gradient wherein random taxa are affected by flow velocity in different ways (described here), and a real-world example wherein a river velocity gradient was sampled at five different points (original data described in Puijalon and Bornette [Reference Puijalon and Bornette2004]).

Simulated Data: River Velocity Gradient

To generate artificial communities simulating the effects of a short gradient, we used a modification of the approach described in Wolda (Reference Wolda1981), wherein a hypothetical fauna is generated using a log-series distribution:

(22)$$M = \alpha \ln \left({1 + \displaystyle{P \over \alpha }} \right).$$

Here, M is the number of species, P is the number of individuals in the population, and α is a diversity coefficient (Fig. 2). We used an α value of 5 with a P of 100,000, corresponding to 50 species (Fig. 2). The relative abundance of each of the 50 species was generated randomly (common and rare species are distributed randomly within the 50 species). To determine the relative abundance of each species, we generated 50 random integers between 1 and 1000; these integers were divided by the sum of the random values, and the resulting string of 50 values was used as a probability vector for randomly selecting the species label a given individual would receive. Next, individual specimens were assigned a randomly generated set of morphological measurements (10 measurements in total)—measurements were generated using R's built-in rnorm function with a mean of 1 and a standard deviation of 3. From this initial community, we generated 9 more communities (representing different points in a hypothetical flow gradient) from which we drew random samples. Each of the 9 new communities was identical to the original, except for the morphological measurements of 10 taxa, selected at random to show one of six morphological trends across the idealized flow gradient (Fig. 2):

1. 1. Minor increase: Two taxa were randomly selected to show a minor increase in measured values along the gradient. Specifically, the mean values for these species increase by 10% between each community.

2. 2. Minor decrease: Two taxa were randomly selected to show a minor decrease in measured values along the gradient. Specifically, the mean values for these species decrease by 10% between each community.

3. 3. Minor fluctuation: Three taxa were randomly selected to show a minor fluctuation in measured values along the gradient. Namely, two were chosen to randomly decrease by 10%, return to the original values, increase by 10%, return to the original values, etc. One taxon was chosen to show the opposite pattern (begin by increasing, etc.).

4. 4. Strong increase: One taxon was chosen to randomly show a strong increase in measured values along the gradient. Specifically, the mean values for this species increase by 20% between each community.

5. 5. Strong decrease: One taxon was chosen to randomly show a strong decrease in measured values along the gradient. Specifically, the mean values for this species decrease by 20% between each community.

6. 6. Subtle decrease: One taxon was chosen to randomly show a very subtle decrease in measured values along the gradient. Specifically, the mean values for this species decrease by 1% between each community.

The remaining 40 taxa were chosen to show no pattern of variation along the gradient. Because the distribution used to generate the random values had a standard deviation of 3, some degree of morphological overlap between communities with close mean values is expected. See R code in the Supplementary Material for more information. The decision to have visible phenotypic effects show up randomly in only a small number of taxa—with a subset showing nonlinear variation—was intended to produce a situation in which randomly choosing a specific taxon that tracks the gradient would be difficult, necessitating the use of an index that works with the entire sample. Furthermore, these underlying populations were designed such that an ordination based on a classic similarity measure (e.g., Jaccard Reference Jaccard1901; Bray and Curtis Reference Bray and Curtis1957) would be unable to reproduce the gradient, as it is contained in the morphology only (Fig. 3).

Figure 2.

Flow diagram illustrating how the hypothetical communities were generated. A, Richness value (here 50 species) is determined based on the size of the underlying population (P), a diversity value (α), and a formula describing a log-series distribution. B, Taxa are randomly generated and assigned random species labels, with some labels being assigned less often (“rare species”). C, A subset of the taxa (n = 10) is randomly chosen to display one of six morphological patterns along the 10 positions of the gradient, while the remaining taxa (n = 40) show no morphological change along the gradient. D, Random samples are pulled from the resulting 10 communities and compared using the proposed indices.

Figure 3.

Ordination plots (principal coordinates analysis) showing the behavior of nine commonly used similarity indices when they are used to compare random samples pulled from the 10 artificial communities described in the text. Each point represents a random sample of 1000 individuals pulled from the corresponding velocity grade. Because the underlying gradient is contained in the morphology of taxa only, none of the common similarity indices produces a plot illustrating the gradient. Incidence-based indices (Jaccard and Sørensen-Dice) produce a plot with multiple samples coinciding on five points; because of their large size, most samples include all the same species (and thus coincide in the plot), with a few samples missing the rarest species in the population by chance alone. The remaining indices are all sensitive to the number of individuals and show random scattering around a central point driven by the random sampling process.