In this document we show how to use the functions described in De Cáceres et al. (2011) by following an example of dietary preferences in pigeons belonging to two populations (Moià and Barcelona) in Catalonia (north east of Spain). We begin by loading the library and the data:
library(indicspecies)
data(pigeons)
The data consists of three items: the resource use matrix for each pigeon population and the matrix of dissimilarities between pairs of resources (seeds).
ls()
#> [1] "dfood" "diet.barcelona" "diet.moia" "groups"
#> [5] "groupskm" "indval" "indvalori" "indvalrest"
#> [9] "indvalspcomb" "p" "pcv" "phi"
#> [13] "prefsign" "prefstat" "sc" "sc2"
#> [17] "wetcomb" "wetkm" "wetland" "wetlandpa"
Before starting any resource niche calculation, we can inspect the matrix of dissimilarities between the \(r=6\) resources, \(\mathbf{D}\):
dfood#> Oats Peas Popcorn Soy Sunflower
#> Peas 0.7029375
#> Popcorn 0.3706196 0.5070959
#> Soy 0.7738801 0.6645980 0.8638219
#> Sunflower 0.7656489 0.9773428 0.8784905 0.8669801
#> Wheat 0.3948828 0.5705974 0.2836482 0.7740283 1.0000000
These can be summarized using a dendrogram. For example:
plot(hclust(dfood, method="average"), h=-1, xlab="",
ylab="Distance", main="", sub="", ylim=c(0.2,1))
Some seeds are quite similar (e.g. popcorn, wheat or oats), whereas sunflower differs substantially from all other resources.
We begin by showing the resource use data for pigeons of Barcelona and Moia, expressed as proportions (i.e., the vector \(\mathbf{f}\) for each population):
= colSums(diet.barcelona)
diet.pop.barcelona round(diet.pop.barcelona/sum(diet.pop.barcelona), dig=3)
#> Oats Peas Popcorn Soy Sunflower Wheat
#> 0.049 0.001 0.028 0.024 0.455 0.442
= colSums(diet.moia)
diet.pop.moia round(diet.pop.moia/sum(diet.pop.moia), dig=3)
#> Oats Peas Popcorn Soy Sunflower Wheat
#> 0.002 0.018 0.091 0.045 0.054 0.791
Whereas pigeons in Moià feed almost exclusively on wheat, those of Barcelona combine wheat and sunflower seeds.
We will determine the resource niche breadth of each of the two populations as De Cáceres et al. (2011):
\(B_{pop} = (1/2)\sum_{j=1}^r\sum_{k=1}^r{f_jf_kd_{jk}}\)
We first conduct our calculations without taking into account the resemblance between resources (which is equal to stating that \(d_{jk}=1\) in all cases):
nichevar(P=diet.barcelona, mode="single")
#> B LC UC
#> Niche 0.2966578 0.2699079 0.3334394
nichevar(P=diet.moia, mode="single")
#> B LC UC
#> Niche 0.1807663 0.1012183 0.2599718
In general, we can say that the niche breadth of the population in Moià is smaller than the niche breadth of the population in Barcelona. If we repeat the same calculations with the matrix of resource resemblance, we realize that the niche breadth of both populations becomes smaller:
= nichevar(P=diet.barcelona, D=dfood,
popvar.barcelona mode="single")
popvar.barcelona#> B LC UC
#> Niche 0.2453822 0.2352728 0.2578303
= nichevar(P=diet.moia, D=dfood, mode="single")
popvar.moia
popvar.moia#> B LC UC
#> Niche 0.0853328 0.02493894 0.1701466
The reason is that the first analysis was assuming that all resources were equally (and maximally) distinct, while the second analysis accounts for the similarity between some resources. Moreover, note that the niche breadth of Moià has decreased more than the niche breadth of Barcelona. This reflects that the resources being used by pigeons of Moià are more similar than the resources used by pigeons of Barcelona.
We can now calculate the niche overlap between the two pigeon populations: \[
O_{12} = \frac{\sum_{j=1}^r\sum_{k=1}^r{f_{1j}f_{2k}(1-d_{jk}^2)}}{\sqrt{\sum_{j=1}^r\sum_{k=1}^r{f_{1j}f_{1k}(1-d_{jk}^2)}\sum_{j=1}^r\sum_{k=1}^r{f_{2j}f_{2k}(1-d_{jk}^2)}}}
\] Using function nicheoverlap()
:
nicheoverlap(P1=diet.barcelona, P2=diet.moia, mode="single")
#> O LC UC
#> Overlap 0.7419319 0.4645122 0.9290476
nicheoverlap(P1=diet.barcelona, P2=diet.moia, mode="single", D = dfood)
#> O LC UC
#> Overlap 0.7912472 0.5411 0.9496495
If we include the resemblance between resources, the degree of overlap increases, for the same reason that we obtained smaller niche breadth statistics when resemblances were included.
In this section, we perform a resource niche analysis at individual level. In particular, we are interested in assessing how much the resource niche of individuals differs from that of their corresponding population. For this, we need to calculate a measure of the degree of individual specialization.
We begin by showing the resource use data for all 23 pigeons in the sample from the Barcelona population, expressed as proportions (i.e., matrix \(\mathbf{F}\) for Barcelona):
round(sweep(diet.barcelona, 1, FUN="/",
rowSums(diet.barcelona)), dig=3)
#> Oats Peas Popcorn Soy Sunflower Wheat
#> 1 0.000 0.000 0.000 0.000 0.020 0.980
#> 2 0.000 0.000 0.000 0.000 1.000 0.000
#> 3 0.000 0.000 0.000 0.000 1.000 0.000
#> 4 0.025 0.000 0.153 0.000 0.000 0.822
#> 5 0.000 0.000 0.000 0.000 0.965 0.035
#> 6 0.000 0.023 0.136 0.000 0.648 0.193
#> 7 0.000 0.000 0.000 0.000 0.795 0.205
#> 8 0.042 0.000 0.000 0.000 0.000 0.958
#> 9 1.000 0.000 0.000 0.000 0.000 0.000
#> 10 0.014 0.000 0.000 0.000 0.014 0.973
#> 11 0.000 0.000 0.000 0.000 0.350 0.650
#> 12 0.000 0.000 0.000 0.000 0.047 0.953
#> 13 0.000 0.000 0.000 0.000 0.437 0.562
#> 14 0.000 0.000 0.000 0.000 0.389 0.611
#> 15 0.000 0.000 0.000 0.549 0.451 0.000
#> 16 0.000 0.000 0.000 0.000 0.913 0.087
#> 17 0.000 0.000 0.113 0.075 0.755 0.057
#> 18 0.000 0.000 0.043 0.000 0.743 0.214
#> 19 0.000 0.000 0.000 0.000 0.667 0.333
#> 20 0.000 0.000 0.000 0.000 0.367 0.633
#> 21 0.000 0.000 0.000 0.000 0.910 0.090
#> 22 0.000 0.000 0.019 0.019 0.786 0.175
#> 23 0.442 0.000 0.000 0.000 0.000 0.558
We see that most individuals in Barcelona feed on either sunflower or wheat, but there are some individuals (like pigeon 9) which prefer oats. Now we display the resource use data for the 19 pigeons representing the population in Moià (i.e., matrix \(\mathbf{F}\) for Moià):
round(sweep(diet.moia, 1, FUN="/",
rowSums(diet.moia)), dig=3)
#> Oats Peas Popcorn Soy Sunflower Wheat
#> 1 0.000 0.000 0.000 0.324 0.000 0.676
#> 2 0.000 0.000 0.000 0.012 0.000 0.988
#> 3 0.000 0.000 0.000 0.000 0.000 1.000
#> 4 0.000 0.000 0.000 0.000 0.000 1.000
#> 5 0.000 0.000 0.000 0.000 0.000 1.000
#> 6 0.000 0.000 0.267 0.000 0.000 0.733
#> 7 0.000 0.000 0.319 0.000 0.000 0.681
#> 8 0.027 0.311 0.041 0.000 0.000 0.622
#> 9 0.000 0.000 0.000 0.000 0.000 1.000
#> 10 0.000 0.000 0.148 0.000 0.000 0.852
#> 11 0.000 0.000 0.200 0.000 0.000 0.800
#> 12 0.000 0.000 0.235 0.000 0.000 0.765
#> 13 0.000 0.000 0.000 0.000 0.000 1.000
#> 14 0.000 0.000 0.119 0.339 0.136 0.407
#> 15 0.000 0.000 0.000 0.000 0.000 1.000
#> 16 0.000 0.000 0.000 0.000 0.000 1.000
#> 17 0.000 0.000 0.430 0.000 0.000 0.570
#> 18 0.000 0.000 0.000 0.000 0.583 0.417
#> 19 0.000 0.000 0.000 0.000 1.000 0.000
Many pigeons from Moià feed on wheat seeds almost exclusively, but some of them have broader preferences.
We begin our resource niche analysis by calculating the niche breadth of each individual in the population, \(B_i\): \[ B_i = (1/2)\sum_{j=1}^r\sum_{k=1}^r{f_{ij}f_{ik}d_{jk}} \] We calculate the values for the individuals of both populations:
<-nichevar(P=diet.barcelona, D=dfood)
indvar.barcelonasummary(indvar.barcelona)
#> B
#> Min. :0.00000
#> 1st Qu.:0.01755
#> Median :0.08153
#> Mean :0.11168
#> 3rd Qu.:0.19890
#> Max. :0.24609
<-nichevar(P=diet.moia, D=dfood)
indvar.moiasummary(indvar.moia)
#> B
#> Min. :0.00000
#> 1st Qu.:0.00000
#> Median :0.01015
#> Mean :0.04031
#> 3rd Qu.:0.01861
#> Max. :0.24306
Most individuals have niche breadths that are smaller than the niche breadth of their corresponding population, although a few individuals in Moià have niche breadths larger than the population value. A niche breadth equal to zero indicates that only one resource is exploited.
We can compare the niche breadth values of the two populations using a non-parametric test.
wilcox.test(indvar.barcelona$B, indvar.moia$B)
#> Warning in wilcox.test.default(indvar.barcelona$B, indvar.moia$B): cannot
#> compute exact p-value with ties
#>
#> Wilcoxon rank sum test with continuity correction
#>
#> data: indvar.barcelona$B and indvar.moia$B
#> W = 328, p-value = 0.005446
#> alternative hypothesis: true location shift is not equal to 0
The Wilcoxon test confirms that the niche breadth of pigeons in Barcelona is generally higher than that of pigeons in Moià, as we saw at the population level.
In order to calculate the degree of individual specialization, Bolnick et al. (2002) defined WIC/TNW, i.e. the ratio between the within individual component (i.e. average niche width) and the total niche width of the population. Similarly we define the following specialization measure, that takes into account the resemblance between resources: \[ S_{pop} = \frac{\sum_{i=n}^n{B_i}/n}{B_{pop}} \] where \(B_i\) is the niche breadth of each individual, and \(B_{pop}\) is the niche breadth of the population. Note that it is possible that \(B_i\) values can be larger than \(B_{pop}\). However, we do not expect the average of \(B_i\) values to be larger than \(B_{pop}\). If we calculate \(S_{pop}\) for the two populations we have:
= mean(indvar.barcelona$B)/popvar.barcelona$B
Spec.barcelona
Spec.barcelona#> [1] 0.4551234
= mean(indvar.moia$B)/popvar.moia$B
Spec.moia
Spec.moia#> [1] 0.4723618
Surprisingly, the degree of specialization in Moià seems slightly higher than that in Barcelona. To see whether this holds statistically, we can calculate the degree of specialization of each individual: \[ S_i = \frac{B_i}{B_{pop}} \] which, in R, is:
= indvar.barcelona$B/popvar.barcelona$B
Spec.ind.barcelona = indvar.moia$B/popvar.moia$B Spec.ind.moia
Finally, we compare this two vectors in a Wilcoxon test:
wilcox.test(Spec.ind.barcelona, Spec.ind.moia)
#> Warning in wilcox.test.default(Spec.ind.barcelona, Spec.ind.moia): cannot
#> compute exact p-value with ties
#>
#> Wilcoxon rank sum test with continuity correction
#>
#> data: Spec.ind.barcelona and Spec.ind.moia
#> W = 277, p-value = 0.1392
#> alternative hypothesis: true location shift is not equal to 0
Which tells us that those differences in individual specialization are not statistically significant.
The idea of this section is to determine how much the niche of each individual overlaps with the niche of other individuals in the population. This can be done by calling function nicheoverlap()
using mode = "pairwise"
:
= nicheoverlap(diet.barcelona,D=dfood, mode="pairwise")
O.barcelona = nicheoverlap(diet.moia,D=dfood, mode="pairwise") O.moia
These calls to nicheoverlap()
return a symmetric square matrix with as many rows and columns as individuals in the resource use data frame. Each cell value in the symmetric matrix is the overlap between two individuals of the population. Using these matrices we can derive the average overlap in each population:
mean(O.barcelona[lower.tri(O.barcelona)])
#> [1] 0.6726374
mean(O.moia[lower.tri(O.moia)])
#> [1] 0.8438816
We can also calculate the average overlap between each individual and the remaining individuals in its population:
= (rowSums(O.barcelona)-1)/(nrow(O.barcelona)-1)
O.barcelona.ind summary(O.barcelona.ind)
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 0.5400 0.6028 0.6558 0.6726 0.7595 0.8055
= (rowSums(O.moia)-1)/(nrow(O.moia)-1)
O.moia.ind summary(O.moia.ind)
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 0.09157 0.89915 0.90649 0.84388 0.90965 0.91134
We substracted one in the numerator and denominator in order to exclude the target individual from the average (the overlap between a resource niche and itself is always one). Apparently, the individuals in Moià have a larger degree of overlap with individuals of their population than individuals in Barcelona. A non-parametric test seems to confirm this difference:
wilcox.test(O.barcelona.ind, O.moia.ind)
#> Warning in wilcox.test.default(O.barcelona.ind, O.moia.ind): cannot compute
#> exact p-value with ties
#>
#> Wilcoxon rank sum test with continuity correction
#>
#> data: O.barcelona.ind and O.moia.ind
#> W = 37, p-value = 4.549e-06
#> alternative hypothesis: true location shift is not equal to 0