Assuming data files are in the working directory
getwd()
## [1] "/Users/lauren/Downloads"
arctic_traits_090424<-as.data.frame(read.csv("arctic_traits_090424_b.csv",
row.names = 1))
View(arctic_traits_090424)
arctic_cpue_090424<-as.data.frame(read.csv("biomass.transformed.matrix.csv",
row.names = 1))
arctic_env_090424<-as.data.frame(read.csv("arctic_env_090424_b.csv",
row.names = 1))
# Six species with no info on age-at-maturity:
arctic_traits_090424[,1:3][rowSums(arctic_traits_090424[,1:3])==0,]
## low.age.at.maturity medium.age.at.maturity
## Anisarchus.medius 0 0
## Enophrys.diceraus 0 0
## Eumesogrammus.praecisus 0 0
## Lycodes.spp. 0 0
## Stichaeus.punctatus 0 0
## Trichocottus.brashnikovi 0 0
## high.age.at.maturity
## Anisarchus.medius 0
## Enophrys.diceraus 0
## Eumesogrammus.praecisus 0
## Lycodes.spp. 0
## Stichaeus.punctatus 0
## Trichocottus.brashnikovi 0
##choose trait taxa from cpue
arctic_traits_090424_b<-as.data.frame(arctic_traits_090424) %>%
filter(rownames(arctic_traits_090424) %in%
colnames(arctic_cpue_090424))
View(arctic_traits_090424_b)
##select stations in arctic_env_090424
arctic_cpue_090424.2 <- as.data.frame(arctic_cpue_090424) %>%
filter(row.names(arctic_cpue_090424) %in%
row.names(arctic_env_090424))
dim(arctic_env_090424)
## [1] 431 7
##create weighted traits matrix
arctic_traits_090424_b<-arctic_traits_090424_b[,-21]
arctic.traits <- as.data.frame(as.matrix(arctic_cpue_090424.2) %*%
as.matrix(arctic_traits_090424_b))
View(arctic.traits)
##create fuzzy coded traits matrix
# 'prep.fuzzy' standardizes modalities within each trait to sum to 1:
arctic.traits.fuzzy<-as.data.frame(prep.fuzzy(arctic.traits,
c(3,6,3,5,3,3)))
# Replace NAs (resulting from traits that have a 0 for all modalities)
arctic.traits.fuzzy[is.na(arctic.traits.fuzzy)] <- 0
arctic.fca <- dudi.fca(arctic.traits.fuzzy, scannf = F, nf=2)
scatter(arctic.fca, csub = 3, clab.moda = 1.5)
summary(arctic.fca)
## Class: fca dudi
## Call: dudi.fca(df = arctic.traits.fuzzy, scannf = F, nf = 2)
##
## Total inertia: 0.1995
##
## Eigenvalues:
## Ax1 Ax2 Ax3 Ax4 Ax5
## 0.07332 0.03414 0.01921 0.01633 0.01144
##
## Projected inertia (%):
## Ax1 Ax2 Ax3 Ax4 Ax5
## 36.758 17.116 9.633 8.187 5.737
##
## Cumulative projected inertia (%):
## Ax1 Ax1:2 Ax1:3 Ax1:4 Ax1:5
## 36.76 53.87 63.51 71.69 77.43
##
## (Only 5 dimensions (out of 17) are shown)
arctic.fca$cr
## RS1 RS2
## FV1 0.05085099 0.008495159
## FV2 0.11368902 0.074642624
## FV3 0.10386336 0.005476967
## FV4 0.06769199 0.064917568
## FV5 0.05066136 0.024185292
## FV6 0.05316375 0.027127278
Cor.all<-cor(arctic_env_090424)
corrplot(Cor.all, method="number", order = "original",
addCoef.col = "black", diag = FALSE, is.corr = TRUE,
type="lower",tl.col = "black", shade.col = "grey50")
Results show moderate correlation between salinity and temperature but do not suggest a need to eliminate any of the variables from the analysis a priori
Conduct hierarchical cluster analysis based on functional traits and Gower distances
fuzzy.dis <- vegdist(arctic.traits.fuzzy, method="gower")
fuzzy.clus <- hclust(fuzzy.dis, method = "ward.D2")
plot(fuzzy.clus, labels=F)
fviz_nbclust(arctic.traits.fuzzy, kmeans,
method = "sil", diss = fuzzy.dis, k.max = 6)
fviz_nbclust(arctic.traits.fuzzy, kmeans,
method = "wss", diss = fuzzy.dis, k.max = 6)
fviz_nbclust(arctic.traits.fuzzy, kmeans,
method = "gap_stat", diss = fuzzy.dis, k.max = 6)
# All methods point towards 3 clusters!
set.seed(1)
fuzzy.clus.kmeans<-kmeans(arctic.traits.fuzzy, centers = 3)
fviz_cluster(fuzzy.clus.kmeans, data = arctic.traits.fuzzy, labelsize=0) +
scale_color_manual(values = c("#0D7EF9", "#98A100", "#FD1C26","orange4")) +
scale_fill_manual(values = c("#0D7EF9", "#98A100", "#FD1C26","orange4"))
There is good separation among the three clusters in just the first two PCs of functional trait composition (CWM matrix), although a PCA based on Euclidean distances is not really an appropriate ordination method here.
## [1] TRUE
Create stacked barplots based on community-weighted means for each modality
Conduct redundancy analysis based on community-weighted means of functional traits. Start with a model including the depth gradient, dynamic environmental variables (temperature, day of ice retreat and salinity) and latitude/longitude to account for any remaining (linear) spatial gradients.
RDA_model.arctic <- rda(arctic.traits.fuzzy ~
Depth + doy_retreat + BottomTemp +
BottomSalinity + Latitude + Longitude,
data = arctic_env)
vif.cca(RDA_model.arctic)
## Depth doy_retreat BottomTemp BottomSalinity Latitude
## 1.185658 1.580370 2.902462 1.998060 1.971497
## Longitude
## 1.577640
RDA_model.arctic
## Call: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp
## + BottomSalinity + Latitude + Longitude, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.04483 0.16746 6
## Unconstrained 0.22289 0.83254 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4 RDA5 RDA6
## 0.03341 0.00877 0.00147 0.00083 0.00027 0.00009
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.07220 0.04100 0.03337 0.02047 0.01305 0.01201 0.00845 0.00634
## (Showing 8 of 17 unconstrained eigenvalues)
RsquareAdj(RDA_model.arctic)
## $r.squared
## [1] 0.1674586
##
## $adj.r.squared
## [1] 0.1556773
plot(RDA_model.arctic)
The full model shows a gradient from warm, low salinity, shallow stations with an early ice retreat to the opposite, and a spatial gradient from SW to NE (Bering Strait to NE Chukchi Sea). This spatial gradient may reflect something like “distance from Bering Strait”(?)
Use permutation-based analysis of variance to assess the statistical significance of the model overall, of RDA axes and of individual terms:
anova(RDA_model.arctic)##significant
## Permutation test for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df Variance F Pr(>F)
## Model 6 0.044832 14.214 0.001 ***
## Residual 424 0.222885
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(RDA_model.arctic, by="axis")
## Permutation test for rda under reduced model
## Forward tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df Variance F Pr(>F)
## RDA1 1 0.033408 63.5531 0.001 ***
## RDA2 1 0.008766 16.6762 0.001 ***
## RDA3 1 0.001468 2.7928 0.250
## RDA4 1 0.000828 1.5758 0.601
## RDA5 1 0.000269 0.5124 0.981
## RDA6 1 0.000091 0.1737 0.992
## Residual 424 0.222885
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(RDA_model.arctic, by="margin")
## Permutation test for rda under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df Variance F Pr(>F)
## Depth 1 0.001999 3.8037 0.004 **
## doy_retreat 1 0.001124 2.1382 0.053 .
## BottomTemp 1 0.002485 4.7277 0.003 **
## BottomSalinity 1 0.000334 0.6349 0.657
## Latitude 1 0.005665 10.7760 0.001 ***
## Longitude 1 0.005650 10.7488 0.001 ***
## Residual 424 0.222885
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
The model is highly significant overall, and the first two RDA axes are significant (intrepretable).
All marginal effects (significance of individual terms after all other terms are included) were significant. The largest variations in trait composition were associated with Temperature and the remaining spatial gradient (latitude/longitude). However, because there are likely other year-to-year differences that are not captured by the covariates, it would be best to use permutations within years for the permutation (rather than across all samples, which is the default). This is similar to considering year as a random effect as residual year-to-year differences are essentially removed by this approach:
anova(RDA_model.arctic, by="margin",
permutations = how(nperm=4999, blocks = arctic_env$Year))
## Permutation test for rda under reduced model
## Marginal effects of terms
## Blocks: arctic_env$Year
## Permutation: free
## Number of permutations: 4999
##
## Model: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df Variance F Pr(>F)
## Depth 1 0.001999 3.8037 0.0200 *
## doy_retreat 1 0.001124 2.1382 0.6354
## BottomTemp 1 0.002485 4.7277 0.0208 *
## BottomSalinity 1 0.000334 0.6349 0.7632
## Latitude 1 0.005665 10.7760 0.0002 ***
## Longitude 1 0.005650 10.7488 0.0002 ***
## Residual 424 0.222885
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
The results are roughly the same, but the day of ice retreat is no longer significant, which is likely because ice retreat differs much more between years than among stations within years, hence it can ‘absorbs’ some of the year to year differences.
plot(doy_retreat ~ factor(Year), data=, arctic_env)
This suggest that a year effect may be appropriate to include in the model, although it needs to be interpreted with caution because of the spatially unbalanced sampling design. Interannual differences in ice retreat are likely confounded with other sources of interannual variability and do not appear to provide the best explanation of such interannual differences when both year and ice retreat are included.
Moreover, ice retreat is also confounded with temperature (r = -0.54), thus its effects are highly uncertain
We conduct a step-wise regression to potentially simplify the model objectively (using forward & backward selection based on adjusted \(R^2\) values).
RDA_model.arctic.2 <- rda( arctic.traits.fuzzy ~ 1,
data = arctic_env)
RDA_mod.step <- ordiR2step(RDA_model.arctic.2,
scope = formula(RDA_model.arctic),
direction="both")
## Step: R2.adj= 0
## Call: arctic.traits.fuzzy ~ 1
##
## R2.adjusted
## <All variables> 0.15567733
## + Latitude 0.10095382
## + BottomTemp 0.08203750
## + Longitude 0.06724984
## + doy_retreat 0.04270619
## + BottomSalinity 0.02466407
## + Depth 0.01524609
## <none> 0.00000000
##
## Df AIC F Pr(>F)
## + Latitude 1 -611.85 49.285 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1009538
## Call: arctic.traits.fuzzy ~ Latitude
##
## R2.adjusted
## <All variables> 0.1556773
## + BottomTemp 0.1246677
## + Longitude 0.1243322
## + Depth 0.1147063
## + BottomSalinity 0.1122931
## + doy_retreat 0.1100879
## <none> 0.1009538
##
## Df AIC F Pr(>F)
## + BottomTemp 1 -622.38 12.622 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1246677
## Call: arctic.traits.fuzzy ~ Latitude + BottomTemp
##
## R2.adjusted
## <All variables> 0.1556773
## + Longitude 0.1480764
## + Depth 0.1319377
## + BottomSalinity 0.1273105
## + doy_retreat 0.1269279
## <none> 0.1246677
##
## Df AIC F Pr(>F)
## + Longitude 1 -633.07 12.76 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1480764
## Call: arctic.traits.fuzzy ~ Latitude + BottomTemp + Longitude
##
## R2.adjusted
## <All variables> 0.1556773
## + Depth 0.1539179
## + doy_retreat 0.1502608
## <none> 0.1480764
## + BottomSalinity 0.1480588
##
## Df AIC F Pr(>F)
## + Depth 1 -635.05 3.9481 0.006 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1539179
## Call: arctic.traits.fuzzy ~ Latitude + BottomTemp + Longitude + Depth
##
## R2.adjusted
## + doy_retreat 0.1564027
## <All variables> 0.1556773
## <none> 0.1539179
## + BottomSalinity 0.1534161
RDA_mod.step
## Call: rda(formula = arctic.traits.fuzzy ~ Latitude + BottomTemp + Longitude
## + Depth, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.04331 0.16179 4
## Unconstrained 0.22440 0.83821 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4
## 0.03330 0.00825 0.00127 0.00050
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.07241 0.04154 0.03340 0.02067 0.01318 0.01235 0.00846 0.00636
## (Showing 8 of 17 unconstrained eigenvalues)
plot(RDA_mod.step)
The model eliminates the ice retreat variable and had a similar \(R^2\) value (16.2%)
To assess how much of the variability can (at most) be explained by depth and dynamic ocean variables, we fit a model without the spatial pattern:
RDA_model.arctic.env <- rda(arctic.traits.fuzzy ~
Depth + doy_retreat + BottomTemp +
BottomSalinity,
data = arctic_env)
RDA_model.arctic.env
## Call: rda(formula = arctic.traits.fuzzy ~ Depth + doy_retreat + BottomTemp
## + BottomSalinity, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.02835 0.10590 4
## Unconstrained 0.23937 0.89410 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4
## 0.023871 0.003214 0.000904 0.000362
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.08286 0.04372 0.03338 0.02078 0.01399 0.01312 0.00846 0.00677
## (Showing 8 of 17 unconstrained eigenvalues)
(R2.1 <- RsquareAdj(RDA_model.arctic.env)$r.squ)
## [1] 0.1058986
About 10.6% of the variability in functional traits can (at most) be accounted for by bathymetry and dynamic ocean variables
To assess how much of the variability can (at most) be explained by the (linear) spatial gradient, we fit a model with ONLY the spatial pattern:
RDA_model.arctic.spatial <- rda(arctic.traits.fuzzy ~ Latitude + Longitude,
data = arctic_env)
RDA_model.arctic.spatial
## Call: rda(formula = arctic.traits.fuzzy ~ Latitude + Longitude, data =
## arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.03438 0.12841 2
## Unconstrained 0.23334 0.87159 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2
## 0.028889 0.005487
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.07614 0.04388 0.03475 0.02103 0.01342 0.01239 0.00858 0.00664
## (Showing 8 of 17 unconstrained eigenvalues)
(R2.2 <- RsquareAdj(RDA_model.arctic.spatial)$r.sq)
## [1] 0.128405
About 12.8% of variability (at most) can be attributed to a linear spatial gradient.
The sum of the \(R^2\) values of the two models above (23.4%) is substantially larger than the variability explained by the full model (15.2%) that includes all variables. The difference is due to multicollinearity; that is, when only the spatial trend is included, it can “absorb” some of the variability otherwise attributed to environmental variables (and vice versa).
Based on these comparisons, we can conclude that environmental covariates may explain up to 10.6% of the variability in trait composition and that environmental covariates, combined with a linear spatial gradient, may account for 16.7% of the variability in trait composition.
We can use conditioning to assess if there is additional interannual variability that is not explained by spatial pattern and oceanographic covariates by conditoning on all variables in the full model and adding year as a factor:
table(arctic_env$Year)
##
## 2007 2008 2009 2010 2012 2015 2017 2019
## 30 15 84 34 37 65 125 41
RDA_model.arctic.year <- rda(arctic.traits.fuzzy ~ factor(Year) +
Condition(Depth + doy_retreat + BottomTemp +
BottomSalinity + Latitude + Longitude),
data = arctic_env)
RDA_model.arctic.year
## Call: rda(formula = arctic.traits.fuzzy ~ factor(Year) + Condition(Depth +
## doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude), data =
## arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Conditional 0.04483 0.16746 6
## Constrained 0.02361 0.08819 7
## Unconstrained 0.19928 0.74435 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4 RDA5 RDA6 RDA7
## 0.010084 0.007932 0.002290 0.001454 0.001354 0.000440 0.000057
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.06661 0.03644 0.02875 0.01636 0.01215 0.01107 0.00784 0.00553
## (Showing 8 of 17 unconstrained eigenvalues)
(R2.3 <- RsquareAdj(RDA_model.arctic.year)$r.sq)
## [1] 0.08819173
anova(RDA_model.arctic.year)
## Permutation test for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ factor(Year) + Condition(Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude), data = arctic_env)
## Df Variance F Pr(>F)
## Model 7 0.02361 7.0581 0.001 ***
## Residual 417 0.19928
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
There are additional significant differences among years that are not accounted for by the full model (including the ice retreat variable) and that account for about 8.8% of the variability in trait composition
If a year effect can reasonably be estimated, we can fit a model with all covariates and use a step-wise approach for model selection to identify the reduced model that maximizes the adjusted \(R^2\) value:
RDA_model.arctic.all <- rda(arctic.traits.fuzzy ~ factor(Year) +
doy_retreat + BottomTemp + BottomSalinity +
Depth + Latitude + Longitude,
data = arctic_env)
vif.cca(RDA_model.arctic.all)
## factor(Year)2008 factor(Year)2009 factor(Year)2010 factor(Year)2012
## 1.591261 3.854139 2.558639 2.301676
## factor(Year)2015 factor(Year)2017 factor(Year)2019 doy_retreat
## 3.117906 4.737821 2.616883 3.729038
## BottomTemp BottomSalinity Depth Latitude
## 3.343987 2.115824 1.247262 2.577273
## Longitude
## 1.729275
RDA_model.arctic.all
## Call: rda(formula = arctic.traits.fuzzy ~ factor(Year) + doy_retreat +
## BottomTemp + BottomSalinity + Depth + Latitude + Longitude, data =
## arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.06844 0.25565 13
## Unconstrained 0.19928 0.74435 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4 RDA5 RDA6 RDA7 RDA8 RDA9 RDA10
## 0.03961 0.01182 0.00864 0.00252 0.00242 0.00132 0.00110 0.00056 0.00024 0.00014
## RDA11 RDA12 RDA13
## 0.00005 0.00001 0.00000
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.06661 0.03644 0.02875 0.01636 0.01215 0.01107 0.00784 0.00553
## (Showing 8 of 17 unconstrained eigenvalues)
RsquareAdj(RDA_model.arctic.all)
## $r.squared
## [1] 0.2556503
##
## $adj.r.squared
## [1] 0.2324452
anova(RDA_model.arctic.all, by="margin",
permutations = how(nperm=4999))
## Permutation test for rda under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 4999
##
## Model: rda(formula = arctic.traits.fuzzy ~ factor(Year) + doy_retreat + BottomTemp + BottomSalinity + Depth + Latitude + Longitude, data = arctic_env)
## Df Variance F Pr(>F)
## factor(Year) 7 0.023610 7.0581 0.0002 ***
## doy_retreat 1 0.001707 3.5710 0.0082 **
## BottomTemp 1 0.002261 4.7313 0.0030 **
## BottomSalinity 1 0.000198 0.4143 0.8590
## Depth 1 0.001544 3.2317 0.0148 *
## Latitude 1 0.005170 10.8180 0.0002 ***
## Longitude 1 0.005716 11.9606 0.0002 ***
## Residual 417 0.199275
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(RDA_model.arctic.all)
Including year in the model results in some of the variables being only marginally significant (ice retreat, depth), but the results seem otherwise consistent with the patterns when a year effect is not included.
step.best <- ordiR2step(RDA_model.arctic.2, scope = formula(RDA_model.arctic.all),
direction="both")
## Step: R2.adj= 0
## Call: arctic.traits.fuzzy ~ 1
##
## R2.adjusted
## <All variables> 0.23244516
## + factor(Year) 0.11001376
## + Latitude 0.10095382
## + BottomTemp 0.08203750
## + Longitude 0.06724984
## + doy_retreat 0.04270619
## + BottomSalinity 0.02466407
## + Depth 0.01524609
## <none> 0.00000000
##
## Df AIC F Pr(>F)
## + factor(Year) 7 -610.29 8.5934 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1100138
## Call: arctic.traits.fuzzy ~ factor(Year)
##
## R2.adjusted
## <All variables> 0.2324452
## + Latitude 0.1851526
## + BottomTemp 0.1724627
## + Longitude 0.1608485
## + doy_retreat 0.1445540
## + BottomSalinity 0.1331800
## + Depth 0.1267104
## <none> 0.1100138
##
## Df AIC F Pr(>F)
## + Latitude 1 -647.33 40.006 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1851526
## Call: arctic.traits.fuzzy ~ factor(Year) + Latitude
##
## R2.adjusted
## <All variables> 0.2324452
## + Longitude 0.2059012
## + BottomTemp 0.2033283
## + Depth 0.1971948
## + BottomSalinity 0.1945589
## + doy_retreat 0.1885190
## <none> 0.1851526
##
## Df AIC F Pr(>F)
## + Longitude 1 -657.47 12.026 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.2059012
## Call: arctic.traits.fuzzy ~ factor(Year) + Latitude + Longitude
##
## R2.adjusted
## <All variables> 0.2324452
## + BottomTemp 0.2244299
## + Depth 0.2172672
## + BottomSalinity 0.2145568
## + doy_retreat 0.2113240
## <none> 0.2059012
##
## Df AIC F Pr(>F)
## + BottomTemp 1 -666.67 11.058 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.2244299
## Call: arctic.traits.fuzzy ~ factor(Year) + Latitude + Longitude + BottomTemp
##
## R2.adjusted
## <All variables> 0.2324452
## + doy_retreat 0.2291436
## + Depth 0.2288895
## <none> 0.2244299
## + BottomSalinity 0.2235796
##
## Df AIC F Pr(>F)
## + doy_retreat 1 -668.32 3.5682 0.004 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.2291436
## Call: arctic.traits.fuzzy ~ factor(Year) + Latitude + Longitude + BottomTemp + doy_retreat
##
## R2.adjusted
## + Depth 0.2335206
## <All variables> 0.2324452
## <none> 0.2291436
## + BottomSalinity 0.2283473
step.best
## Call: rda(formula = arctic.traits.fuzzy ~ factor(Year) + Latitude +
## Longitude + BottomTemp + doy_retreat, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 0.26772 1.00000
## Constrained 0.06662 0.24886 11
## Unconstrained 0.20109 0.75114 17
##
## Inertia is variance
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4 RDA5 RDA6 RDA7 RDA8 RDA9 RDA10
## 0.03928 0.01116 0.00845 0.00251 0.00230 0.00126 0.00097 0.00047 0.00017 0.00005
## RDA11
## 0.00001
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 0.06708 0.03706 0.02876 0.01637 0.01215 0.01125 0.00797 0.00553
## (Showing 8 of 17 unconstrained eigenvalues)
anova(step.best, by="margin")
## Permutation test for rda under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ factor(Year) + Latitude + Longitude + BottomTemp + doy_retreat, data = arctic_env)
## Df Variance F Pr(>F)
## factor(Year) 7 0.024281 7.2276 0.001 ***
## Latitude 1 0.005191 10.8163 0.001 ***
## Longitude 1 0.006486 13.5147 0.001 ***
## BottomTemp 1 0.005140 10.7089 0.001 ***
## doy_retreat 1 0.001713 3.5682 0.007 **
## Residual 419 0.201092
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(step.best, by="axis")
## Permutation test for rda under reduced model
## Forward tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = arctic.traits.fuzzy ~ factor(Year) + Latitude + Longitude + BottomTemp + doy_retreat, data = arctic_env)
## Df Variance F Pr(>F)
## RDA1 1 0.039280 81.8457 0.001 ***
## RDA2 1 0.011161 23.2559 0.001 ***
## RDA3 1 0.008450 17.6070 0.001 ***
## RDA4 1 0.002507 5.2235 0.096 .
## RDA5 1 0.002300 4.7914 0.099 .
## RDA6 1 0.001255 2.6152 0.564
## RDA7 1 0.000968 2.0171 0.718
## RDA8 1 0.000474 0.9885 0.978
## RDA9 1 0.000168 0.3508 0.999
## RDA10 1 0.000051 0.1062 1.000
## RDA11 1 0.000010 0.0198
## Residual 419 0.201092
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(step.best)
(R2.fun.best <- RsquareAdj(step.best)$r.sq)
## [1] 0.2488631
The resulting ‘best’ model eliminates day of ice retreat and depth and explains approximately 24.9% of overall variability. It clearly suggests that bottom temperature is a major driver of trait composition. In addition, there is a strong remaining spatial gradient, as well as considerable year-to-year variability.
Plot ordination of functional traits with individual modalities and environmental variables overlaid.
## RDA1 RDA2
## 14.79 4.41
## cumulative contributions of most influential species:
##
## $Intermediate_Southern
## fusiform...normal very.high.longevity medium.high.fecundity
## 0.0700453 0.1368361 0.1954519
## juvenile medium.body.size very.high.fecundity
## 0.2529954 0.3096603 0.3643692
## small.body.size flat medium.longevity
## 0.4189860 0.4712934 0.5230296
## high.age.at.maturity medium.fecundity sub.adult
## 0.5737627 0.6238831 0.6705561
## eel.like
## 0.7155096
##
## $Intermediate_Northern
## very.high.longevity fusiform...normal juvenile
## 0.08500226 0.15082343 0.21580073
## eel.like medium.high.fecundity small.body.size
## 0.27892946 0.33938269 0.39969286
## medium.longevity medium.fecundity medium.body.size
## 0.45861152 0.51526972 0.56902686
## high.longevity sub.adult adult1
## 0.61623060 0.66079559 0.70041729
##
## $Southern_Northern
## very.high.longevity juvenile small.body.size
## 0.1052709 0.1825539 0.2533913
## medium.longevity medium.body.size high.age.at.maturity
## 0.3233804 0.3740326 0.4235775
## sub.adult high.longevity very.high.fecundity
## 0.4703233 0.5156984 0.5590295
## eel.like medium.fecundity flat
## 0.6023426 0.6423036 0.6816705
## large.body.size
## 0.7206724
##
## Contrast: Intermediate_Southern
##
## average sd ratio ava avb cumsum p
## fusiform...normal 0.014822 0.011310 1.310500 0.469200 0.315400 0.070 0.01
## very.high.longevity 0.014133 0.011613 1.217000 0.357800 0.507100 0.137 1.00
## medium.high.fecundity 0.012403 0.010149 1.222200 0.382500 0.275000 0.196 0.02
## juvenile 0.012177 0.008748 1.392000 0.617100 0.727000 0.253 1.00
## medium.body.size 0.011991 0.009730 1.232300 0.565100 0.640200 0.310 0.54
## very.high.fecundity 0.011577 0.008540 1.355600 0.046300 0.177700 0.364 0.01
## small.body.size 0.011557 0.008369 1.381000 0.261200 0.164500 0.419 1.00
## flat 0.011069 0.007287 1.519000 0.037600 0.162600 0.471 0.01
## medium.longevity 0.010948 0.007887 1.388100 0.322400 0.222400 0.523 1.00
## high.age.at.maturity 0.010735 0.007192 1.492600 0.108000 0.222700 0.574 0.01
## medium.fecundity 0.010606 0.008038 1.319400 0.281200 0.336000 0.624 0.90
## sub.adult 0.009876 0.007167 1.378000 0.269000 0.176500 0.671 0.62
## eel.like 0.009512 0.007772 1.223900 0.376200 0.416200 0.716 1.00
## medium.age.at.maturity 0.009450 0.007265 1.300700 0.699200 0.630500 0.760 0.05
## large.body.size 0.008828 0.006742 1.309500 0.173600 0.195300 0.802 0.19
## high.longevity 0.007561 0.006524 1.159000 0.319800 0.270500 0.838 1.00
## low.medium.fecundity 0.007211 0.005204 1.385700 0.116300 0.054000 0.872 0.97
## low.age.at.maturity 0.006942 0.005767 1.203800 0.192800 0.146800 0.904 1.00
## high.fecundity 0.006550 0.005080 1.289300 0.164800 0.152300 0.935 0.23
## adult1 0.006065 0.004457 1.360700 0.113900 0.096600 0.964 1.00
## deep.and.or.short 0.005005 0.004131 1.211400 0.105000 0.093700 0.988 0.98
## compressed 0.001531 0.001805 0.848000 0.012100 0.012100 0.995 0.01
## low.fecundity 0.001056 0.001947 0.542700 0.009000 0.005000 1.000 1.00
##
## fusiform...normal **
## very.high.longevity
## medium.high.fecundity *
## juvenile
## medium.body.size
## very.high.fecundity **
## small.body.size
## flat **
## medium.longevity
## high.age.at.maturity **
## medium.fecundity
## sub.adult
## eel.like
## medium.age.at.maturity *
## large.body.size
## high.longevity
## low.medium.fecundity
## low.age.at.maturity
## high.fecundity
## adult1
## deep.and.or.short
## compressed **
## low.fecundity
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Intermediate_Northern
##
## average sd ratio ava avb cumsum p
## very.high.longevity 0.018713 0.011361 1.647100 0.357800 0.156200 0.085 0.01
## fusiform...normal 0.014490 0.010393 1.394300 0.469200 0.316500 0.151 0.01
## juvenile 0.014305 0.010877 1.315100 0.617100 0.472800 0.216 0.04
## eel.like 0.013898 0.009773 1.422000 0.376200 0.527800 0.279 0.01
## medium.high.fecundity 0.013309 0.009902 1.344100 0.382500 0.260200 0.339 0.01
## small.body.size 0.013277 0.010715 1.239100 0.261200 0.396200 0.400 0.13
## medium.longevity 0.012971 0.010257 1.264600 0.322400 0.446400 0.459 0.05
## medium.fecundity 0.012473 0.009135 1.365500 0.281200 0.395800 0.515 0.01
## medium.body.size 0.011835 0.009921 1.192900 0.565100 0.505800 0.569 0.67
## high.longevity 0.010392 0.008591 1.209600 0.319800 0.397400 0.616 0.01
## sub.adult 0.009811 0.007666 1.279800 0.269000 0.318700 0.661 0.71
## adult1 0.008723 0.006506 1.340800 0.113900 0.208500 0.700 0.01
## medium.age.at.maturity 0.008454 0.006792 1.244700 0.699200 0.706800 0.739 0.98
## low.age.at.maturity 0.008355 0.006731 1.241300 0.192800 0.221900 0.777 0.02
## large.body.size 0.008139 0.006352 1.281400 0.173600 0.098000 0.814 0.95
## low.medium.fecundity 0.007641 0.006106 1.251400 0.116300 0.166500 0.848 0.41
## high.age.at.maturity 0.007276 0.005070 1.435000 0.108000 0.071300 0.881 1.00
## high.fecundity 0.006654 0.005175 1.285900 0.164800 0.101600 0.912 0.17
## deep.and.or.short 0.005724 0.004947 1.157200 0.105000 0.107900 0.938 0.10
## very.high.fecundity 0.004821 0.004508 1.069300 0.046300 0.041400 0.960 1.00
## flat 0.004719 0.004820 0.979100 0.037600 0.043700 0.981 1.00
## low.fecundity 0.002966 0.003363 0.881900 0.009000 0.034500 0.995 0.01
## compressed 0.001204 0.001879 0.641000 0.012100 0.004100 1.000 0.90
##
## very.high.longevity **
## fusiform...normal **
## juvenile *
## eel.like **
## medium.high.fecundity **
## small.body.size
## medium.longevity *
## medium.fecundity **
## medium.body.size
## high.longevity **
## sub.adult
## adult1 **
## medium.age.at.maturity
## low.age.at.maturity *
## large.body.size
## low.medium.fecundity
## high.age.at.maturity
## high.fecundity
## deep.and.or.short .
## very.high.fecundity
## flat
## low.fecundity **
## compressed
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Southern_Northern
##
## average sd ratio ava avb cumsum p
## very.high.longevity 0.029450 0.013988 2.105500 0.507100 0.156200 0.105 0.01
## juvenile 0.021621 0.012562 1.721100 0.727000 0.472800 0.183 0.01
## small.body.size 0.019817 0.012428 1.594600 0.164500 0.396200 0.253 0.01
## medium.longevity 0.019580 0.012083 1.620400 0.222400 0.446400 0.323 0.01
## medium.body.size 0.014170 0.011617 1.219800 0.640200 0.505800 0.374 0.01
## high.age.at.maturity 0.013861 0.008229 1.684300 0.222700 0.071300 0.424 0.01
## sub.adult 0.013078 0.008926 1.465100 0.176500 0.318700 0.470 0.01
## high.longevity 0.012694 0.010043 1.264000 0.270500 0.397400 0.516 0.01
## very.high.fecundity 0.012122 0.008612 1.407500 0.177700 0.041400 0.559 0.01
## eel.like 0.012117 0.008401 1.442400 0.416200 0.527800 0.602 0.01
## medium.fecundity 0.011179 0.008436 1.325200 0.336000 0.395800 0.642 0.23
## flat 0.011013 0.007369 1.494600 0.162600 0.043700 0.682 0.01
## large.body.size 0.010911 0.007737 1.410300 0.195300 0.098000 0.721 0.01
## medium.age.at.maturity 0.010608 0.008449 1.255500 0.630500 0.706800 0.759 0.01
## low.medium.fecundity 0.010485 0.007006 1.496600 0.054000 0.166500 0.796 0.01
## medium.high.fecundity 0.010298 0.007879 1.307100 0.275000 0.260200 0.833 1.00
## adult1 0.010125 0.007197 1.406900 0.096600 0.208500 0.869 0.01
## fusiform...normal 0.009932 0.007716 1.287100 0.315400 0.316500 0.905 1.00
## low.age.at.maturity 0.009488 0.007607 1.247300 0.146800 0.221900 0.939 0.01
## high.fecundity 0.007183 0.005197 1.382100 0.152300 0.101600 0.964 0.01
## deep.and.or.short 0.005842 0.004870 1.199600 0.093700 0.107900 0.985 0.02
## low.fecundity 0.002961 0.003555 0.832800 0.005000 0.034500 0.996 0.01
## compressed 0.001223 0.002088 0.586000 0.012100 0.004100 1.000 0.82
##
## very.high.longevity **
## juvenile **
## small.body.size **
## medium.longevity **
## medium.body.size **
## high.age.at.maturity **
## sub.adult **
## high.longevity **
## very.high.fecundity **
## eel.like **
## medium.fecundity
## flat **
## large.body.size **
## medium.age.at.maturity **
## low.medium.fecundity **
## medium.high.fecundity
## adult1 **
## fusiform...normal
## low.age.at.maturity **
## high.fecundity **
## deep.and.or.short *
## low.fecundity **
## compressed
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 99
arctic_cpue <- as.data.frame(arctic_cpue_090424.2) %>%
filter(row.names(arctic_cpue_090424.2) %in%
row.names(arctic_env_090424))
# Bray-Curtis distances:
tax.dis <- vegdist(arctic_cpue, method="bray")
tax.clus <- hclust(tax.dis, method = "ward.D2")
plot(tax.clus, labels=F)
There are two obvious clusters, but it is unclear if the 2nd cluster should be grouped into more smaller clusters and, if so, into how many. There appears to be a fairly obvious cutoff at 5 clusters (at a height of about 3.5)
fviz_nbclust(arctic_cpue, kmeans,
method = "silh", diss = tax.dis, k.max = 8) # 2 (or 4)
fviz_nbclust(arctic_cpue, kmeans,
method = "wss", diss = tax.dis, k.max = 8) # 2 or 5?
fviz_nbclust(arctic_cpue, kmeans,
method = "gap", diss = tax.dis, k.max = 8) # 5
There is some ambiguity but results suggest either 2 or 5 clusters, consistent with the dendrogram. Nevertheless, for consistency with the functional clusters and for ease of comparison, we use three clusters:
set.seed(1)
tax.clus.kmeans.1 <- kmeans(arctic_cpue, centers = 4,
iter.max = 100)
# Or, if using 5 clusters (messier!)
# tax.clus.kmeans<-kmeans(arctic_cpue, centers = 5)
# Visualize clusters using a PCA of the taxonomic data (which is unlikely to represent data very well)
fviz_cluster(tax.clus.kmeans.1, data = arctic_cpue, labelsize=0) +
scale_color_manual(values = c("#FD1C26", "#98A100", "#0D7EF9","orange")) +
scale_fill_manual(values = c("#FD1C26", "#98A100", "#0D7EF9","orange"))
Compared to the functional clusters, the PCA of the cluster results accounts for a limited proportion of the variability in species abundances and the clusters overlap substantially in the space of the first two PCs. This could be a result of the PCA not being an appropriate ordination method for the species abundance data (many zeros!).
## [1] TRUE
There was a comment from someone (Caitlin?): “I added distance to shore in addition to latitude and longitude as our spatial component we are specifically thinking about the Alaska Coastal Current compared to Anadyr so distance to shore in tandem with temp and sal should capture that?”
Franz: This makes sense, but Distance from shore was not in the data frame!
Also, I made a few changes below. I used the reduced data set (23 taxa) and used cca instead of rda, as cca allows for a Gaussian (dome-shaped) response of taxa to environmental gradients.
CCA_model.arctic.tax <- cca(arctic_cpue ~
factor(Year) + Depth + doy_retreat + BottomTemp +
BottomSalinity + Latitude + Longitude,
data = arctic_env)
vif.cca(CCA_model.arctic.tax)
## factor(Year)2008 factor(Year)2009 factor(Year)2010 factor(Year)2012
## 1.604783 3.314902 1.849274 2.135118
## factor(Year)2015 factor(Year)2017 factor(Year)2019 Depth
## 2.484255 3.941739 2.137603 1.373861
## doy_retreat BottomTemp BottomSalinity Latitude
## 3.527987 3.586086 2.558954 2.452021
## Longitude
## 1.807384
CCA_model.arctic.tax
## Call: cca(formula = arctic_cpue ~ factor(Year) + Depth + doy_retreat +
## BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 2.6272 1.0000
## Constrained 0.4912 0.1870 13
## Unconstrained 2.1360 0.8130 32
##
## Inertia is scaled Chi-square
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## CCA1 CCA2 CCA3 CCA4 CCA5 CCA6 CCA7 CCA8 CCA9 CCA10
## 0.17185 0.08012 0.07115 0.04144 0.03229 0.02698 0.01985 0.01472 0.01289 0.00748
## CCA11 CCA12 CCA13
## 0.00668 0.00316 0.00257
##
## Eigenvalues for unconstrained axes:
## CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
## 0.16433 0.12777 0.11988 0.10552 0.10031 0.09188 0.08794 0.08743
## (Showing 8 of 32 unconstrained eigenvalues)
(R2.4 <- RsquareAdj(CCA_model.arctic.tax)$r.sq)
## [1] 0.1869642
plot(CCA_model.arctic.tax)
Similar to the functional traits CWM-RDA analysis, the results show a gradient from warm, low salinity, shallow stations with an early ice retreat to the opposite, and a spatial gradient from SW to NE (Bering Strait to NE Chukchi Sea).
Assess significance of model overall, CCA axes and individual terms (marginal effects)
anova(CCA_model.arctic.tax)
## Permutation test for cca under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: cca(formula = arctic_cpue ~ factor(Year) + Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df ChiSquare F Pr(>F)
## Model 13 0.49118 7.3764 0.001 ***
## Residual 417 2.13597
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(CCA_model.arctic.tax, by="axis") ##CCA1-7 significant
## Permutation test for cca under reduced model
## Forward tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: cca(formula = arctic_cpue ~ factor(Year) + Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df ChiSquare F Pr(>F)
## CCA1 1 0.17185 33.5494 0.001 ***
## CCA2 1 0.08012 15.6411 0.001 ***
## CCA3 1 0.07115 13.8904 0.001 ***
## CCA4 1 0.04144 8.0912 0.001 ***
## CCA5 1 0.03229 6.3043 0.001 ***
## CCA6 1 0.02698 5.2672 0.001 ***
## CCA7 1 0.01985 3.8759 0.011 *
## CCA8 1 0.01472 2.8747 0.306
## CCA9 1 0.01289 2.5166 0.674
## CCA10 1 0.00748 1.4600 1.000
## CCA11 1 0.00668 1.3033
## CCA12 1 0.00316 0.6175
## CCA13 1 0.00257 0.5010
## Residual 417 2.13597
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(CCA_model.arctic.tax, by="margin")
## Permutation test for cca under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: cca(formula = arctic_cpue ~ factor(Year) + Depth + doy_retreat + BottomTemp + BottomSalinity + Latitude + Longitude, data = arctic_env)
## Df ChiSquare F Pr(>F)
## factor(Year) 7 0.20697 5.7724 0.001 ***
## Depth 1 0.03094 6.0400 0.001 ***
## doy_retreat 1 0.01194 2.3310 0.001 ***
## BottomTemp 1 0.02124 4.1463 0.001 ***
## BottomSalinity 1 0.00897 1.7513 0.048 *
## Latitude 1 0.03760 7.3404 0.001 ***
## Longitude 1 0.02461 4.8036 0.001 ***
## Residual 417 2.13597
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
The overall model is significant and accounts for approximately 18.7% of the variability in species composition.
The first 5 CCA axes capture significant proportions of overall variability.
All individual terms except salinity are significant.
CCA_model.arctic.2.tax <- cca( arctic_cpue ~ 1 ,
data = arctic_env)
step.res.tax<-ordiR2step(CCA_model.arctic.2.tax,
scope = formula(CCA_model.arctic.tax),
direction="both")
## Step: R2.adj= 0
## Call: arctic_cpue ~ 1
##
## R2.adjusted
## <All variables> 0.16058896
## + factor(Year) 0.07685339
## + BottomTemp 0.04456610
## + Latitude 0.03966945
## + BottomSalinity 0.02988888
## + Depth 0.02575235
## + Longitude 0.01878999
## + doy_retreat 0.01744261
## <none> 0.00000000
##
## Df AIC F Pr(>F)
## + factor(Year) 7 1995.7 6.1574 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.07679814
## Call: arctic_cpue ~ factor(Year)
##
## R2.adjusted
## <All variables> 0.16058896
## + BottomTemp 0.11985292
## + Latitude 0.11144451
## + BottomSalinity 0.10666028
## + Depth 0.10356078
## + doy_retreat 0.09172164
## + Longitude 0.09130296
## <none> 0.07679814
##
## Df AIC F Pr(>F)
## + BottomTemp 1 1976.1 21.699 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1198166
## Call: arctic_cpue ~ factor(Year) + BottomTemp
##
## R2.adjusted
## <All variables> 0.1605890
## + Latitude 0.1334319
## + Longitude 0.1319710
## + Depth 0.1315062
## + BottomSalinity 0.1260115
## + doy_retreat 0.1246509
## <none> 0.1198166
##
## Df AIC F Pr(>F)
## + Latitude 1 1970.3 7.711 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.133516
## Call: arctic_cpue ~ factor(Year) + BottomTemp + Latitude
##
## R2.adjusted
## <All variables> 0.1605890
## + Depth 0.1468811
## + Longitude 0.1445690
## + BottomSalinity 0.1398196
## + doy_retreat 0.1374759
## <none> 0.1335160
##
## Df AIC F Pr(>F)
## + Depth 1 1964.5 7.6592 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1469695
## Call: arctic_cpue ~ factor(Year) + BottomTemp + Latitude + Depth
##
## R2.adjusted
## <All variables> 0.1605890
## + Longitude 0.1563189
## + BottomSalinity 0.1503971
## + doy_retreat 0.1499596
## <none> 0.1469695
##
## Df AIC F Pr(>F)
## + Longitude 1 1960.7 5.7389 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1563866
## Call: arctic_cpue ~ factor(Year) + BottomTemp + Latitude + Depth + Longitude
##
## R2.adjusted
## <All variables> 0.1605890
## + doy_retreat 0.1592537
## + BottomSalinity 0.1579646
## <none> 0.1563866
##
## Df AIC F Pr(>F)
## + doy_retreat 1 1960.2 2.4277 0.002 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1591602
## Call: arctic_cpue ~ factor(Year) + BottomTemp + Latitude + Depth + Longitude + doy_retreat
##
## R2.adjusted
## <All variables> 0.1605890
## + BottomSalinity 0.1605691
## <none> 0.1591602
##
## Df AIC F Pr(>F)
## + BottomSalinity 1 1960.3 1.7513 0.036 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: R2.adj= 0.1605687
## Call: arctic_cpue ~ factor(Year) + BottomTemp + Latitude + Depth + Longitude + doy_retreat + BottomSalinity
##
step.res.tax
## Call: cca(formula = arctic_cpue ~ factor(Year) + BottomTemp + Latitude +
## Depth + Longitude + doy_retreat + BottomSalinity, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 2.6272 1.0000
## Constrained 0.4912 0.1870 13
## Unconstrained 2.1360 0.8130 32
##
## Inertia is scaled Chi-square
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## CCA1 CCA2 CCA3 CCA4 CCA5 CCA6 CCA7 CCA8 CCA9 CCA10
## 0.17185 0.08012 0.07115 0.04144 0.03229 0.02698 0.01985 0.01472 0.01289 0.00748
## CCA11 CCA12 CCA13
## 0.00668 0.00316 0.00257
##
## Eigenvalues for unconstrained axes:
## CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
## 0.16433 0.12777 0.11988 0.10552 0.10031 0.09188 0.08794 0.08743
## (Showing 8 of 32 unconstrained eigenvalues)
anova(step.res.tax, by="margin")
## Permutation test for cca under reduced model
## Marginal effects of terms
## Permutation: free
## Number of permutations: 999
##
## Model: cca(formula = arctic_cpue ~ factor(Year) + BottomTemp + Latitude + Depth + Longitude + doy_retreat + BottomSalinity, data = arctic_env)
## Df ChiSquare F Pr(>F)
## factor(Year) 7 0.20697 5.7724 0.001 ***
## BottomTemp 1 0.02124 4.1463 0.001 ***
## Latitude 1 0.03760 7.3404 0.001 ***
## Depth 1 0.03094 6.0400 0.001 ***
## Longitude 1 0.02461 4.8036 0.001 ***
## doy_retreat 1 0.01194 2.3310 0.003 **
## BottomSalinity 1 0.00897 1.7513 0.044 *
## Residual 417 2.13597
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(R2.best <- RsquareAdj(step.res.tax)$r.sq)
## [1] 0.1869642
plot(step.res.tax)
The ‘best’ model with the highest adjusted \(R^2\) eliminates salinity only but retains all other variables and all of the remaining variables are highly significant. It explains slightly less of the overall variability in species composition than the full model (18.7%)
We fit a model without the spatial pattern to assess how much of the variability can (at most) be explained by depth and dynamic ocean variables:
CCA_model.arctic.env.tax <- cca(arctic_cpue ~
Depth + doy_retreat + BottomTemp +
BottomSalinity,
data = arctic_env)
CCA_model.arctic.env.tax
## Call: cca(formula = arctic_cpue ~ Depth + doy_retreat + BottomTemp +
## BottomSalinity, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 2.62715 1.00000
## Constrained 0.20527 0.07813 4
## Unconstrained 2.42188 0.92187 32
##
## Inertia is scaled Chi-square
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## CCA1 CCA2 CCA3 CCA4
## 0.13664 0.04110 0.01552 0.01201
##
## Eigenvalues for unconstrained axes:
## CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
## 0.18807 0.14467 0.14086 0.12835 0.11414 0.11089 0.10317 0.10123
## (Showing 8 of 32 unconstrained eigenvalues)
(R2.5 <- RsquareAdj(CCA_model.arctic.env.tax)$r.sq)
## [1] 0.07813255
About 7.8% of variability (at most) can be accounted for by bathymetry and dynamic ocean variables
Next, we fit a model with ONLY the spatial pattern to assess how much of the variability can (at most) be explained by the (linear) spatial gradient:
CCA_model.arctic.spatial.tax <- cca(arctic_cpue ~ Latitude + Longitude,
data = arctic_env)
CCA_model.arctic.spatial.tax
## Call: cca(formula = arctic_cpue ~ Latitude + Longitude, data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 2.6271 1.0000
## Constrained 0.1497 0.0570 2
## Unconstrained 2.4774 0.9430 32
##
## Inertia is scaled Chi-square
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## CCA1 CCA2
## 0.11098 0.03877
##
## Eigenvalues for unconstrained axes:
## CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
## 0.18696 0.17954 0.14929 0.12910 0.12519 0.11219 0.10270 0.10017
## (Showing 8 of 32 unconstrained eigenvalues)
(R2.6 <- RsquareAdj(CCA_model.arctic.spatial.tax)$r.sq)
## [1] 0.05700268
Only about 5.7% of variability (at most) can be attributed to a linear spatial gradient.
Given that the overall model (including year) explains about 18.7% of the variability, there appears to be considerable year-to-year variability not explained by any of the environmental or spatial gradients.
The additional interannual variability that is not explained by spatial pattern and oceanographic covariates can be quantified by conditoning on all other variables in the best model:
step.res.tax2 <- cca(arctic_cpue ~ factor(Year) + Condition(BottomTemp +
Latitude + Depth + doy_retreat + Longitude),
data = arctic_env)
step.res.tax2
## Call: cca(formula = arctic_cpue ~ factor(Year) + Condition(BottomTemp +
## Latitude + Depth + doy_retreat + Longitude), data = arctic_env)
##
## -- Model Summary --
##
## Inertia Proportion Rank
## Total 2.62715 1.00000
## Conditional 0.27208 0.10356 5
## Constrained 0.21014 0.07999 7
## Unconstrained 2.14494 0.81645 32
##
## Inertia is scaled Chi-square
##
## -- Eigenvalues --
##
## Eigenvalues for constrained axes:
## CCA1 CCA2 CCA3 CCA4 CCA5 CCA6 CCA7
## 0.07257 0.04729 0.03567 0.02491 0.01740 0.00642 0.00588
##
## Eigenvalues for unconstrained axes:
## CA1 CA2 CA3 CA4 CA5 CA6 CA7 CA8
## 0.16490 0.12912 0.12029 0.10552 0.10071 0.09236 0.08915 0.08776
## (Showing 8 of 32 unconstrained eigenvalues)
anova(step.res.tax2)
## Permutation test for cca under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: cca(formula = arctic_cpue ~ factor(Year) + Condition(BottomTemp + Latitude + Depth + doy_retreat + Longitude), data = arctic_env)
## Df ChiSquare F Pr(>F)
## Model 7 0.21014 5.8502 0.001 ***
## Residual 418 2.14494
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(R2.7 <- RsquareAdj(step.res.tax2)$r.sq)
## [1] 0.07998691
There are additional significant differences among years that are not accounted for by other variables in the best model that account for about 8% of the variability in taxonomic composition
Plot CCA ordination with species and environmental variables overlaid. (this seems to use different scaling for environmental vectors from default CCA plot?)
## [1] "BottomTemp" "Latitude" "Depth" "Longitude" "doy_retreat"
## CCA1 CCA2
## 6.54 3.05
##
## Contrast: Nearshore_Intermediate
##
## average sd ratio ava avb cumsum
## Myoxocephalus.scorpius 0.05948 0.03392 1.75320 6.48000 1.85100 0.090
## Gymnocanthus.tricuspis 0.05697 0.03755 1.51700 6.54300 2.38200 0.176
## Artediellus.scaber 0.04323 0.03514 1.23010 3.74000 0.70700 0.242
## Boreogadus.saida 0.04253 0.03551 1.19740 4.49300 3.32000 0.306
## Lumpenus.fabricii.sp. 0.03781 0.03302 1.14520 6.01600 4.40300 0.363
## Stichaeus.punctatus 0.03190 0.03290 0.96980 2.63800 0.76000 0.411
## Liparis.spp. 0.03190 0.03394 0.93970 1.44600 2.19600 0.460
## Aspidophoroides.olrikii 0.03013 0.02921 1.03130 2.26200 1.62200 0.505
## Hippoglossoides.robustus 0.02870 0.02924 0.98160 1.38900 2.05600 0.549
## Eleginus.gracilis 0.02798 0.03339 0.83810 2.25100 0.63200 0.591
## Lycodes.spp. 0.02333 0.02689 0.86750 1.26700 1.42300 0.626
## Nautichthys.pribilovius 0.02331 0.02660 0.87630 1.85100 0.44200 0.662
## Triglops.pingelii 0.02133 0.02691 0.79250 1.71900 0.38100 0.694
## Liparis.tunicatus 0.02082 0.02800 0.74330 1.11600 1.04600 0.725
## Podothecus.veternus 0.01867 0.02769 0.67430 1.49200 0.14400 0.754
## Anisarchus.medius 0.01649 0.02513 0.65610 0.59800 0.99700 0.779
## Limanda.aspera 0.01542 0.02659 0.58000 1.18300 0.28700 0.802
## Gymnelus.hemifasciatus 0.01408 0.02520 0.55880 0.91900 0.47300 0.823
## Liparis.gibbus 0.01369 0.02341 0.58460 0.86100 0.49100 0.844
## Trichocottus.brashnikovi 0.01234 0.02188 0.56400 0.96500 0.20800 0.863
## Ammodytes.hexapterus 0.01185 0.02463 0.48130 0.87000 0.19100 0.881
## Gymnelus.viridis 0.01045 0.02065 0.50590 0.66700 0.32300 0.896
## Icelus.spatula 0.01026 0.02202 0.46600 0.31300 0.63700 0.912
## Gadus.chalcogrammus 0.01025 0.02133 0.48050 0.33300 0.59100 0.927
## Pleuronectidae 0.00769 0.01992 0.38580 0.35800 0.34600 0.939
## Gadus.macrocephalus 0.00691 0.01835 0.37680 0.29100 0.34400 0.950
## Eumesogrammus.praecisus 0.00671 0.01783 0.37630 0.46600 0.15600 0.960
## Podothecus.accipenserinus 0.00550 0.01604 0.34260 0.37600 0.07700 0.968
## Aspidophoroides.monopterygius 0.00538 0.01406 0.38280 0.46700 0.07400 0.976
## Hemilepidotus.papilio 0.00488 0.01577 0.30960 0.39400 0.04200 0.984
## Limanda.proboscidea 0.00462 0.01442 0.32020 0.29500 0.11600 0.991
## Gymnelus.sp. 0.00390 0.01322 0.29530 0.17600 0.17900 0.996
## Leptoclinus.maculatus 0.00239 0.01046 0.22850 0.14800 0.08900 1.000
## p
## Myoxocephalus.scorpius 0.01 **
## Gymnocanthus.tricuspis 0.01 **
## Artediellus.scaber 0.01 **
## Boreogadus.saida 0.01 **
## Lumpenus.fabricii.sp. 1.00
## Stichaeus.punctatus 0.01 **
## Liparis.spp. 0.02 *
## Aspidophoroides.olrikii 1.00
## Hippoglossoides.robustus 1.00
## Eleginus.gracilis 0.01 **
## Lycodes.spp. 1.00
## Nautichthys.pribilovius 0.01 **
## Triglops.pingelii 0.01 **
## Liparis.tunicatus 0.34
## Podothecus.veternus 0.01 **
## Anisarchus.medius 1.00
## Limanda.aspera 0.01 **
## Gymnelus.hemifasciatus 1.00
## Liparis.gibbus 0.38
## Trichocottus.brashnikovi 0.01 **
## Ammodytes.hexapterus 0.05 *
## Gymnelus.viridis 0.41
## Icelus.spatula 0.98
## Gadus.chalcogrammus 0.33
## Pleuronectidae 0.01 **
## Gadus.macrocephalus 0.01 **
## Eumesogrammus.praecisus 0.02 *
## Podothecus.accipenserinus 0.01 **
## Aspidophoroides.monopterygius 0.92
## Hemilepidotus.papilio 0.05 *
## Limanda.proboscidea 0.09 .
## Gymnelus.sp. 0.06 .
## Leptoclinus.maculatus 1.00
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Nearshore_Offshore
##
## average sd ratio ava avb cumsum
## Anisarchus.medius 0.04419 0.02419 1.82690 0.59800 4.97000 0.083
## Hippoglossoides.robustus 0.04029 0.02479 1.62560 1.38900 5.15600 0.158
## Lycodes.spp. 0.03334 0.02236 1.49150 1.26700 4.15900 0.220
## Artediellus.scaber 0.03215 0.02636 1.21990 3.74000 1.78900 0.280
## Myoxocephalus.scorpius 0.03023 0.02487 1.21560 6.48000 4.13800 0.337
## Boreogadus.saida 0.02972 0.02491 1.19310 4.49300 5.14800 0.393
## Aspidophoroides.olrikii 0.02755 0.02180 1.26420 2.26200 3.75300 0.444
## Stichaeus.punctatus 0.02449 0.02523 0.97080 2.63800 0.43000 0.490
## Lumpenus.fabricii.sp. 0.02264 0.02098 1.07920 6.01600 6.51300 0.532
## Eleginus.gracilis 0.02263 0.02524 0.89670 2.25100 0.96300 0.575
## Liparis.spp. 0.02127 0.02486 0.85570 1.44600 1.56800 0.614
## Gymnocanthus.tricuspis 0.01951 0.01582 1.23380 6.54300 5.74800 0.651
## Liparis.tunicatus 0.01829 0.02227 0.82130 1.11600 1.43900 0.685
## Nautichthys.pribilovius 0.01773 0.01994 0.88930 1.85100 0.43300 0.718
## Triglops.pingelii 0.01760 0.02027 0.86810 1.71900 0.78600 0.751
## Podothecus.veternus 0.01711 0.02144 0.79810 1.49200 0.84800 0.783
## Gymnelus.hemifasciatus 0.01434 0.02048 0.70010 0.91900 1.01700 0.810
## Liparis.gibbus 0.01185 0.01792 0.66110 0.86100 0.72700 0.832
## Limanda.aspera 0.01174 0.02032 0.57790 1.18300 0.23200 0.854
## Aspidophoroides.monopterygius 0.01029 0.01658 0.62030 0.46700 0.90600 0.873
## Ammodytes.hexapterus 0.00959 0.01916 0.50070 0.87000 0.24700 0.891
## Trichocottus.brashnikovi 0.00933 0.01689 0.55230 0.96500 0.14000 0.909
## Gymnelus.viridis 0.00656 0.01482 0.44240 0.66700 0.10900 0.921
## Gadus.chalcogrammus 0.00566 0.01326 0.42680 0.33300 0.36100 0.932
## Leptoclinus.maculatus 0.00528 0.01314 0.40170 0.14800 0.48200 0.942
## Hemilepidotus.papilio 0.00523 0.01367 0.38270 0.39400 0.21800 0.951
## Icelus.spatula 0.00482 0.01251 0.38510 0.31300 0.27200 0.960
## Eumesogrammus.praecisus 0.00480 0.01356 0.35390 0.46600 0.09900 0.969
## Gadus.macrocephalus 0.00387 0.01272 0.30450 0.29100 0.14900 0.977
## Podothecus.accipenserinus 0.00363 0.01159 0.31350 0.37600 0.00000 0.983
## Limanda.proboscidea 0.00355 0.01064 0.33340 0.29500 0.10700 0.990
## Pleuronectidae 0.00354 0.01237 0.28630 0.35800 0.03700 0.997
## Gymnelus.sp. 0.00186 0.00786 0.23710 0.17600 0.03000 1.000
## p
## Anisarchus.medius 0.02 *
## Hippoglossoides.robustus 0.01 **
## Lycodes.spp. 1.00
## Artediellus.scaber 0.80
## Myoxocephalus.scorpius 1.00
## Boreogadus.saida 1.00
## Aspidophoroides.olrikii 1.00
## Stichaeus.punctatus 0.01 **
## Lumpenus.fabricii.sp. 1.00
## Eleginus.gracilis 0.01 **
## Liparis.spp. 1.00
## Gymnocanthus.tricuspis 1.00
## Liparis.tunicatus 0.91
## Nautichthys.pribilovius 0.01 **
## Triglops.pingelii 0.01 **
## Podothecus.veternus 0.01 **
## Gymnelus.hemifasciatus 1.00
## Liparis.gibbus 0.86
## Limanda.aspera 0.01 **
## Aspidophoroides.monopterygius 0.01 **
## Ammodytes.hexapterus 0.46
## Trichocottus.brashnikovi 0.01 **
## Gymnelus.viridis 1.00
## Gadus.chalcogrammus 1.00
## Leptoclinus.maculatus 0.84
## Hemilepidotus.papilio 0.01 **
## Icelus.spatula 1.00
## Eumesogrammus.praecisus 0.32
## Gadus.macrocephalus 0.50
## Podothecus.accipenserinus 0.03 *
## Limanda.proboscidea 0.31
## Pleuronectidae 0.53
## Gymnelus.sp. 0.77
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Nearshore_Northern
##
## average sd ratio ava avb cumsum
## Gymnocanthus.tricuspis 0.07042 0.03349 2.10290 6.54300 0.79100 0.098
## Myoxocephalus.scorpius 0.06259 0.03246 1.92820 6.48000 1.34200 0.185
## Lumpenus.fabricii.sp. 0.05546 0.03568 1.55410 6.01600 1.80100 0.262
## Anisarchus.medius 0.05354 0.03062 1.74890 0.59800 4.90000 0.336
## Lycodes.spp. 0.04496 0.03004 1.49680 1.26700 4.49100 0.399
## Artediellus.scaber 0.03940 0.03278 1.20180 3.74000 2.01700 0.454
## Boreogadus.saida 0.03396 0.02888 1.17610 4.49300 5.39500 0.501
## Aspidophoroides.olrikii 0.03028 0.02816 1.07530 2.26200 2.17600 0.543
## Stichaeus.punctatus 0.03007 0.03265 0.92090 2.63800 0.03200 0.585
## Eleginus.gracilis 0.02535 0.03288 0.77110 2.25100 0.13400 0.620
## Hippoglossoides.robustus 0.02420 0.02729 0.88680 1.38900 1.53400 0.654
## Liparis.spp. 0.02410 0.03126 0.77100 1.44600 1.16000 0.687
## Nautichthys.pribilovius 0.02170 0.02544 0.85280 1.85100 0.22400 0.718
## Gymnelus.hemifasciatus 0.02124 0.02717 0.78200 0.91900 1.49500 0.747
## Triglops.pingelii 0.02024 0.02567 0.78840 1.71900 0.33200 0.775
## Podothecus.veternus 0.01738 0.02625 0.66210 1.49200 0.05700 0.799
## Liparis.tunicatus 0.01715 0.02636 0.65070 1.11600 0.57900 0.823
## Icelus.spatula 0.01527 0.02355 0.64830 0.31300 1.24100 0.845
## Ammodytes.hexapterus 0.01386 0.02546 0.54420 0.87000 0.48000 0.864
## Limanda.aspera 0.01352 0.02486 0.54380 1.18300 0.10100 0.882
## Liparis.gibbus 0.01351 0.02224 0.60770 0.86100 0.56600 0.901
## Gymnelus.viridis 0.01275 0.02156 0.59140 0.66700 0.69800 0.919
## Trichocottus.brashnikovi 0.01076 0.02060 0.52240 0.96500 0.02900 0.934
## Gadus.chalcogrammus 0.00733 0.01717 0.42700 0.33300 0.37200 0.944
## Eumesogrammus.praecisus 0.00608 0.01668 0.36470 0.46600 0.14500 0.953
## Leptoclinus.maculatus 0.00552 0.01466 0.37620 0.14800 0.41500 0.960
## Aspidophoroides.monopterygius 0.00549 0.01454 0.37800 0.46700 0.08400 0.968
## Podothecus.accipenserinus 0.00473 0.01475 0.32100 0.37600 0.02500 0.975
## Hemilepidotus.papilio 0.00428 0.01458 0.29350 0.39400 0.00000 0.981
## Limanda.proboscidea 0.00411 0.01358 0.30270 0.29500 0.07300 0.986
## Pleuronectidae 0.00395 0.01466 0.26970 0.35800 0.00000 0.992
## Gadus.macrocephalus 0.00320 0.01352 0.23690 0.29100 0.00000 0.996
## Gymnelus.sp. 0.00274 0.01112 0.24630 0.17600 0.07200 1.000
## p
## Gymnocanthus.tricuspis 0.01 **
## Myoxocephalus.scorpius 0.01 **
## Lumpenus.fabricii.sp. 0.01 **
## Anisarchus.medius 0.01 **
## Lycodes.spp. 0.01 **
## Artediellus.scaber 0.01 **
## Boreogadus.saida 0.97
## Aspidophoroides.olrikii 1.00
## Stichaeus.punctatus 0.01 **
## Eleginus.gracilis 0.01 **
## Hippoglossoides.robustus 1.00
## Liparis.spp. 1.00
## Nautichthys.pribilovius 0.01 **
## Gymnelus.hemifasciatus 0.32
## Triglops.pingelii 0.01 **
## Podothecus.veternus 0.01 **
## Liparis.tunicatus 0.96
## Icelus.spatula 0.33
## Ammodytes.hexapterus 0.01 **
## Limanda.aspera 0.01 **
## Liparis.gibbus 0.49
## Gymnelus.viridis 0.03 *
## Trichocottus.brashnikovi 0.01 **
## Gadus.chalcogrammus 0.95
## Eumesogrammus.praecisus 0.05 *
## Leptoclinus.maculatus 0.84
## Aspidophoroides.monopterygius 0.96
## Podothecus.accipenserinus 0.01 **
## Hemilepidotus.papilio 0.04 *
## Limanda.proboscidea 0.09 .
## Pleuronectidae 0.53
## Gadus.macrocephalus 0.77
## Gymnelus.sp. 0.26
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Intermediate_Offshore
##
## average sd ratio ava avb cumsum
## Anisarchus.medius 0.05374 0.03214 1.67210 0.99700 4.97000 0.091
## Hippoglossoides.robustus 0.04645 0.03269 1.42100 2.05600 5.15600 0.169
## Gymnocanthus.tricuspis 0.04640 0.03233 1.43510 2.38200 5.74800 0.248
## Boreogadus.saida 0.04261 0.03454 1.23360 3.32000 5.14800 0.320
## Myoxocephalus.scorpius 0.04212 0.03185 1.32250 1.85100 4.13800 0.391
## Lycodes.spp. 0.04195 0.02909 1.44210 1.42300 4.15900 0.462
## Aspidophoroides.olrikii 0.03769 0.02806 1.34330 1.62200 3.75300 0.525
## Lumpenus.fabricii.sp. 0.03505 0.03156 1.11050 4.40300 6.51300 0.585
## Liparis.spp. 0.03119 0.03154 0.98880 2.19600 1.56800 0.637
## Artediellus.scaber 0.02440 0.03010 0.81060 0.70700 1.78900 0.678
## Liparis.tunicatus 0.02286 0.02766 0.82630 1.04600 1.43900 0.717
## Gymnelus.hemifasciatus 0.01594 0.02665 0.59830 0.47300 1.01700 0.744
## Eleginus.gracilis 0.01594 0.02452 0.65020 0.63200 0.96300 0.771
## Liparis.gibbus 0.01254 0.02158 0.58120 0.49100 0.72700 0.792
## Stichaeus.punctatus 0.01243 0.02253 0.55150 0.76000 0.43000 0.813
## Triglops.pingelii 0.01179 0.02044 0.57700 0.38100 0.78600 0.833
## Podothecus.veternus 0.01110 0.02194 0.50590 0.14400 0.84800 0.852
## Aspidophoroides.monopterygius 0.01057 0.01993 0.53040 0.07400 0.90600 0.870
## Gadus.chalcogrammus 0.01055 0.02136 0.49390 0.59100 0.36100 0.887
## Icelus.spatula 0.00995 0.02127 0.46770 0.63700 0.27200 0.904
## Nautichthys.pribilovius 0.00902 0.01855 0.48640 0.44200 0.43300 0.919
## Leptoclinus.maculatus 0.00651 0.01695 0.38400 0.08900 0.48200 0.930
## Limanda.aspera 0.00581 0.01696 0.34270 0.28700 0.23200 0.940
## Gadus.macrocephalus 0.00559 0.01587 0.35210 0.34400 0.14900 0.950
## Gymnelus.viridis 0.00502 0.01491 0.33630 0.32300 0.10900 0.958
## Ammodytes.hexapterus 0.00471 0.01464 0.32180 0.19100 0.24700 0.966
## Pleuronectidae 0.00452 0.01546 0.29210 0.34600 0.03700 0.974
## Trichocottus.brashnikovi 0.00356 0.01180 0.30190 0.20800 0.14000 0.980
## Hemilepidotus.papilio 0.00321 0.01249 0.25740 0.04200 0.21800 0.985
## Eumesogrammus.praecisus 0.00306 0.01156 0.26460 0.15600 0.09900 0.990
## Limanda.proboscidea 0.00253 0.01162 0.21800 0.11600 0.10700 0.995
## Gymnelus.sp. 0.00226 0.01023 0.22110 0.17900 0.03000 0.999
## Podothecus.accipenserinus 0.00092 0.00631 0.14520 0.07700 0.00000 1.000
## p
## Anisarchus.medius 0.01 **
## Hippoglossoides.robustus 0.01 **
## Gymnocanthus.tricuspis 0.01 **
## Boreogadus.saida 0.01 **
## Myoxocephalus.scorpius 0.22
## Lycodes.spp. 0.01 **
## Aspidophoroides.olrikii 0.01 **
## Lumpenus.fabricii.sp. 1.00
## Liparis.spp. 0.10 .
## Artediellus.scaber 1.00
## Liparis.tunicatus 0.03 *
## Gymnelus.hemifasciatus 1.00
## Eleginus.gracilis 0.68
## Liparis.gibbus 0.71
## Stichaeus.punctatus 0.99
## Triglops.pingelii 0.95
## Podothecus.veternus 0.48
## Aspidophoroides.monopterygius 0.01 **
## Gadus.chalcogrammus 0.26
## Icelus.spatula 1.00
## Nautichthys.pribilovius 0.99
## Leptoclinus.maculatus 0.57
## Limanda.aspera 0.98
## Gadus.macrocephalus 0.07 .
## Gymnelus.viridis 1.00
## Ammodytes.hexapterus 1.00
## Pleuronectidae 0.22
## Trichocottus.brashnikovi 0.98
## Hemilepidotus.papilio 0.33
## Eumesogrammus.praecisus 0.92
## Limanda.proboscidea 0.73
## Gymnelus.sp. 0.63
## Podothecus.accipenserinus 0.97
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Intermediate_Northern
##
## average sd ratio ava avb cumsum
## Anisarchus.medius 0.07063 0.04353 1.62270 0.99700 4.90000 0.101
## Lycodes.spp. 0.06166 0.04320 1.42750 1.42300 4.49100 0.190
## Lumpenus.fabricii.sp. 0.06053 0.04706 1.28630 4.40300 1.80100 0.276
## Boreogadus.saida 0.05560 0.04715 1.17910 3.32000 5.39500 0.356
## Hippoglossoides.robustus 0.04021 0.04164 0.96550 2.05600 1.53400 0.414
## Liparis.spp. 0.04015 0.04346 0.92390 2.19600 1.16000 0.472
## Aspidophoroides.olrikii 0.03985 0.03936 1.01230 1.62200 2.17600 0.529
## Gymnocanthus.tricuspis 0.03976 0.03975 1.00030 2.38200 0.79100 0.586
## Myoxocephalus.scorpius 0.03761 0.04085 0.92080 1.85100 1.34200 0.640
## Artediellus.scaber 0.03544 0.04144 0.85530 0.70700 2.01700 0.691
## Gymnelus.hemifasciatus 0.02713 0.03913 0.69330 0.47300 1.49500 0.729
## Icelus.spatula 0.02444 0.03537 0.69080 0.63700 1.24100 0.764
## Liparis.tunicatus 0.02195 0.03341 0.65700 1.04600 0.57900 0.796
## Liparis.gibbus 0.01482 0.02852 0.51980 0.49100 0.56600 0.817
## Gadus.chalcogrammus 0.01481 0.03001 0.49360 0.59100 0.37200 0.838
## Gymnelus.viridis 0.01403 0.02751 0.51000 0.32300 0.69800 0.859
## Stichaeus.punctatus 0.01241 0.02784 0.44590 0.76000 0.03200 0.876
## Eleginus.gracilis 0.01136 0.02636 0.43100 0.63200 0.13400 0.893
## Triglops.pingelii 0.01008 0.02321 0.43450 0.38100 0.33200 0.907
## Ammodytes.hexapterus 0.00977 0.02581 0.37870 0.19100 0.48000 0.921
## Nautichthys.pribilovius 0.00942 0.02310 0.40770 0.44200 0.22400 0.935
## Leptoclinus.maculatus 0.00718 0.02023 0.35470 0.08900 0.41500 0.945
## Limanda.aspera 0.00568 0.01987 0.28610 0.28700 0.10100 0.953
## Pleuronectidae 0.00544 0.01986 0.27370 0.34600 0.00000 0.961
## Gadus.macrocephalus 0.00527 0.01794 0.29380 0.34400 0.00000 0.968
## Eumesogrammus.praecisus 0.00438 0.01568 0.27930 0.15600 0.14500 0.975
## Gymnelus.sp. 0.00352 0.01522 0.23160 0.17900 0.07200 0.980
## Podothecus.veternus 0.00329 0.01580 0.20840 0.14400 0.05700 0.985
## Trichocottus.brashnikovi 0.00314 0.01275 0.24660 0.20800 0.02900 0.989
## Limanda.proboscidea 0.00288 0.01593 0.18060 0.11600 0.07300 0.993
## Aspidophoroides.monopterygius 0.00266 0.01382 0.19270 0.07400 0.08400 0.997
## Podothecus.accipenserinus 0.00150 0.00904 0.16540 0.07700 0.02500 0.999
## Hemilepidotus.papilio 0.00065 0.00644 0.10170 0.04200 0.00000 1.000
## p
## Anisarchus.medius 0.01 **
## Lycodes.spp. 0.01 **
## Lumpenus.fabricii.sp. 0.01 **
## Boreogadus.saida 0.01 **
## Hippoglossoides.robustus 0.01 **
## Liparis.spp. 0.01 **
## Aspidophoroides.olrikii 0.01 **
## Gymnocanthus.tricuspis 1.00
## Myoxocephalus.scorpius 1.00
## Artediellus.scaber 0.01 **
## Gymnelus.hemifasciatus 0.01 **
## Icelus.spatula 0.01 **
## Liparis.tunicatus 0.10 .
## Liparis.gibbus 0.13
## Gadus.chalcogrammus 0.01 **
## Gymnelus.viridis 0.01 **
## Stichaeus.punctatus 0.96
## Eleginus.gracilis 1.00
## Triglops.pingelii 1.00
## Ammodytes.hexapterus 0.44
## Nautichthys.pribilovius 0.99
## Leptoclinus.maculatus 0.23
## Limanda.aspera 0.99
## Pleuronectidae 0.01 **
## Gadus.macrocephalus 0.06 .
## Eumesogrammus.praecisus 0.40
## Gymnelus.sp. 0.06 .
## Podothecus.veternus 1.00
## Trichocottus.brashnikovi 1.00
## Limanda.proboscidea 0.69
## Aspidophoroides.monopterygius 1.00
## Podothecus.accipenserinus 0.88
## Hemilepidotus.papilio 1.00
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Contrast: Offshore_Northern
##
## average sd ratio ava avb cumsum
## Gymnocanthus.tricuspis 0.06051 0.02664 2.27130 5.74800 0.79100 0.109
## Lumpenus.fabricii.sp. 0.05940 0.03456 1.71870 6.51300 1.80100 0.215
## Hippoglossoides.robustus 0.04909 0.03155 1.55580 5.15600 1.53400 0.304
## Myoxocephalus.scorpius 0.04240 0.03156 1.34330 4.13800 1.34200 0.380
## Aspidophoroides.olrikii 0.03418 0.02735 1.24960 3.75300 2.17600 0.441
## Artediellus.scaber 0.02952 0.03000 0.98400 1.78900 2.01700 0.494
## Boreogadus.saida 0.02920 0.02750 1.06150 5.14800 5.39500 0.547
## Anisarchus.medius 0.02807 0.02631 1.06670 4.97000 4.90000 0.597
## Lycodes.spp. 0.02738 0.02636 1.03860 4.15900 4.49100 0.646
## Liparis.spp. 0.02394 0.02834 0.84480 1.56800 1.16000 0.689
## Gymnelus.hemifasciatus 0.02225 0.02756 0.80740 1.01700 1.49500 0.729
## Liparis.tunicatus 0.01969 0.02633 0.74760 1.43900 0.57900 0.764
## Icelus.spatula 0.01503 0.02306 0.65180 0.27200 1.24100 0.791
## Liparis.gibbus 0.01242 0.02050 0.60590 0.72700 0.56600 0.814
## Eleginus.gracilis 0.01192 0.02205 0.54080 0.96300 0.13400 0.835
## Triglops.pingelii 0.01094 0.01920 0.56960 0.78600 0.33200 0.855
## Aspidophoroides.monopterygius 0.01041 0.01970 0.52850 0.90600 0.08400 0.874
## Podothecus.veternus 0.00987 0.02017 0.48960 0.84800 0.05700 0.891
## Leptoclinus.maculatus 0.00895 0.01841 0.48620 0.48200 0.41500 0.907
## Gymnelus.viridis 0.00822 0.01805 0.45560 0.10900 0.69800 0.922
## Ammodytes.hexapterus 0.00774 0.01899 0.40770 0.24700 0.48000 0.936
## Gadus.chalcogrammus 0.00765 0.01732 0.44160 0.36100 0.37200 0.950
## Nautichthys.pribilovius 0.00673 0.01573 0.42770 0.43300 0.22400 0.962
## Stichaeus.punctatus 0.00503 0.01519 0.33140 0.43000 0.03200 0.971
## Limanda.aspera 0.00362 0.01313 0.27600 0.23200 0.10100 0.978
## Hemilepidotus.papilio 0.00265 0.01124 0.23600 0.21800 0.00000 0.982
## Eumesogrammus.praecisus 0.00252 0.01004 0.25150 0.09900 0.14500 0.987
## Limanda.proboscidea 0.00209 0.01073 0.19460 0.10700 0.07300 0.991
## Gadus.macrocephalus 0.00177 0.00931 0.19010 0.14900 0.00000 0.994
## Trichocottus.brashnikovi 0.00173 0.00827 0.20930 0.14000 0.02900 0.997
## Gymnelus.sp. 0.00109 0.00728 0.14910 0.03000 0.07200 0.999
## Pleuronectidae 0.00047 0.00493 0.09450 0.03700 0.00000 1.000
## Podothecus.accipenserinus 0.00023 0.00272 0.08470 0.00000 0.02500 1.000
## p
## Gymnocanthus.tricuspis 0.01 **
## Lumpenus.fabricii.sp. 0.01 **
## Hippoglossoides.robustus 0.01 **
## Myoxocephalus.scorpius 0.17
## Aspidophoroides.olrikii 0.17
## Artediellus.scaber 0.99
## Boreogadus.saida 1.00
## Anisarchus.medius 1.00
## Lycodes.spp. 1.00
## Liparis.spp. 1.00
## Gymnelus.hemifasciatus 0.13
## Liparis.tunicatus 0.56
## Icelus.spatula 0.48
## Liparis.gibbus 0.83
## Eleginus.gracilis 1.00
## Triglops.pingelii 1.00
## Aspidophoroides.monopterygius 0.01 **
## Podothecus.veternus 0.88
## Leptoclinus.maculatus 0.01 **
## Gymnelus.viridis 0.93
## Ammodytes.hexapterus 0.91
## Gadus.chalcogrammus 0.94
## Nautichthys.pribilovius 1.00
## Stichaeus.punctatus 1.00
## Limanda.aspera 1.00
## Hemilepidotus.papilio 0.69
## Eumesogrammus.praecisus 0.98
## Limanda.proboscidea 0.92
## Gadus.macrocephalus 0.99
## Trichocottus.brashnikovi 1.00
## Gymnelus.sp. 0.98
## Pleuronectidae 1.00
## Podothecus.accipenserinus 1.00
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 99
Look at diversity metrics (RaoQ and Simpson’s diversity) compared to
clusters
fd.chukchi<-dbFD(arctic_traits_090424_b,
arctic_cpue_090424.2,
w.abun = TRUE,corr = "none",
calc.FRic = TRUE)
## FEVe: Could not be calculated for communities with <3 functionally singular species.
## FRic: To respect s > t, FRic could not be calculated for communities with <3 functionally singular species.
## FRic: Dimensionality reduction was required. The last 20 PCoA axes (out of 22 in total) were removed.
## FRic: Quality of the reduced-space representation = 0.3694546
## FDiv: Could not be calculated for communities with <3 functionally singular species.
View(as.data.frame(fd.chukchi$FDis))
arctic_env_fd<-arctic_env_tax
summary(arctic_env_fd)
## Year Latitude Longitude doy_retreat
## Min. :2007 Min. :65.65 Min. :-169.6 Min. :119.0
## 1st Qu.:2009 1st Qu.:69.50 1st Qu.:-166.9 1st Qu.:143.0
## Median :2015 Median :70.90 Median :-165.0 Median :161.0
## Mean :2013 Mean :70.35 Mean :-164.8 Mean :163.6
## 3rd Qu.:2017 3rd Qu.:71.39 3rd Qu.:-163.0 3rd Qu.:183.0
## Max. :2019 Max. :73.00 Max. :-154.4 Max. :233.0
## Depth BottomTemp BottomSalinity cluster
## Min. :15.00 Min. :-1.760 Min. :27.15 Southern :133
## 1st Qu.:39.90 1st Qu.: 0.328 1st Qu.:32.06 Intermediate:170
## Median :43.00 Median : 2.140 Median :32.28 Northern :128
## Mean :42.87 Mean : 2.600 Mean :32.18
## 3rd Qu.:47.00 3rd Qu.: 4.439 3rd Qu.:32.50
## Max. :90.00 Max. :10.660 Max. :34.56
## tax.cluster
## Nearshore : 86
## Intermediate: 95
## Offshore :111
## Northern :139
##
##
arctic_env_fd$raoQ<-as.numeric(fd.chukchi$RaoQ)
simpson<-as.data.frame(diversity(arctic_cpue, "simpson"))
names(simpson)<-"simpson"
arctic_env_fd$simpson<-simpson$simpson
library(glmm)
## Warning: package 'glmm' was built under R version 4.3.3
## Loading required package: trust
## Loading required package: mvtnorm
## Loading required package: Matrix
##
## Attaching package: 'Matrix'
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
## Loading required package: parallel
## Loading required package: doParallel
## Loading required package: foreach
##
## Attaching package: 'foreach'
## The following objects are masked from 'package:purrr':
##
## accumulate, when
## Loading required package: iterators
library(mgcv)
## Loading required package: nlme
##
## Attaching package: 'nlme'
## The following object is masked from 'package:dplyr':
##
## collapse
## This is mgcv 1.8-42. For overview type 'help("mgcv-package")'.
lm.fd<-lm(raoQ~ cluster,
data=arctic_env_fd)
lm.fd.1<-lm(raoQ~1,
data=arctic_env_fd)
anova(lm.fd, lm.fd.1)
## Analysis of Variance Table
##
## Model 1: raoQ ~ cluster
## Model 2: raoQ ~ 1
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 428 3084.3
## 2 430 3181.0 -2 -96.714 6.7103 0.001351 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
fd.aov<-aov(lm.fd)
summary(lm.fd)
##
## Call:
## lm(formula = raoQ ~ cluster, data = arctic_env_fd)
##
## Residuals:
## Min 1Q Median 3Q Max
## -12.9821 -1.5045 0.3644 1.8041 6.3515
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.3232 0.2328 87.309 < 2e-16 ***
## clusterIntermediate -0.2099 0.3108 -0.675 0.499822
## clusterNorthern -1.1367 0.3324 -3.420 0.000687 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.684 on 428 degrees of freedom
## Multiple R-squared: 0.0304, Adjusted R-squared: 0.02587
## F-statistic: 6.71 on 2 and 428 DF, p-value: 0.001351
TukeyHSD(fd.aov)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = lm.fd)
##
## $cluster
## diff lwr upr p adj
## Intermediate-Southern -0.2098722 -0.9407545 0.5210101 0.7779705
## Northern-Southern -1.1366887 -1.9184350 -0.3549424 0.0019820
## Northern-Intermediate -0.9268165 -1.6656645 -0.1879685 0.0093702
library(rcompanion)
## Warning: package 'rcompanion' was built under R version 4.3.3
Sum = groupwiseMean(raoQ ~ cluster,
data = arctic_env_fd,
conf = 0.95,
digits = 3)
clus.cols.fuzzy <- c("#FD1C26","#98A100", "#0D7EF9")
fuzzy.cluster.CI<-
ggplot(data = Sum,
aes(x=cluster,y= Mean, color=cluster)) +
geom_errorbar(aes(ymin=Trad.lower,
ymax=Trad.upper,
width=0.15))+
geom_point(shape = 15,
size = 4) +
scale_color_manual(values=clus.cols.fuzzy)+
theme_bw()+
labs(x="Cluster", y="Rao's Quadratic Entropy")
fuzzy.cluster.CI
##taxonomy
lm.tax<-lm(simpson~ tax.cluster,
data=arctic_env_fd)
summary(lm.tax)
##
## Call:
## lm(formula = simpson ~ tax.cluster, data = arctic_env_fd)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.32198 -0.01993 0.00883 0.03188 0.11637
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.884585 0.006214 142.344 < 2e-16 ***
## tax.clusterIntermediate -0.075286 0.008578 -8.777 < 2e-16 ***
## tax.clusterOffshore 0.005333 0.008279 0.644 0.52
## tax.clusterNorthern -0.062609 0.007906 -7.919 2.1e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.05763 on 427 degrees of freedom
## Multiple R-squared: 0.278, Adjusted R-squared: 0.2729
## F-statistic: 54.8 on 3 and 427 DF, p-value: < 2.2e-16
tax.aov<-aov(lm.tax)
summary(tax.aov)
## Df Sum Sq Mean Sq F value Pr(>F)
## tax.cluster 3 0.546 0.18199 54.8 <2e-16 ***
## Residuals 427 1.418 0.00332
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(tax.aov)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = lm.tax)
##
## $tax.cluster
## diff lwr upr p adj
## Intermediate-Nearshore -0.075285514 -0.097409213 -0.05316182 0.0000000
## Offshore-Nearshore 0.005332905 -0.016019749 0.02668556 0.9175544
## Northern-Nearshore -0.062609293 -0.083001494 -0.04221709 0.0000000
## Offshore-Intermediate 0.080618419 0.059843479 0.10139336 0.0000000
## Northern-Intermediate 0.012676221 -0.007110243 0.03246269 0.3504801
## Northern-Offshore -0.067942198 -0.086862610 -0.04902179 0.0000000
Sum.simpson = groupwiseMean(simpson ~ tax.cluster,
data = arctic_env_fd,
conf = 0.95,
digits = 3)
clus.cols.tax <- c("#FD1C26","#98A100", "#0D7EF9","orange")
Sum.simpson$tax.cluster <- factor(Sum.simpson$tax.cluster ,
levels = c("Nearshore",
"Intermediate",
"Offshore",
"Northern"))
tax.cluster.CI<-
ggplot(data = Sum.simpson,
aes(x=tax.cluster,y= Mean, color=tax.cluster)) +
geom_errorbar(aes(ymin=Trad.lower,
ymax=Trad.upper,
width=0.15))+
geom_point(shape = 15,
size = 4) +
scale_color_manual(values=clus.cols)+
theme_bw()+
labs(x="Cluster", y="Simpson Diversity Index")
tax.cluster.CI
Diversity metrics compared to environmental gradients for functional
diversity
arctic.cwm.rda.sites$raoQ<-as.numeric(arctic_env_fd$raoQ)
lm.rda.0<-lm(raoQ~ 1,
data=arctic.cwm.rda.sites)
lm.rda.1<-lm(raoQ~ RDA1 + I(RDA1^2) + RDA2 + I(RDA2^2),
data=arctic.cwm.rda.sites)
anova(lm.rda.0, lm.rda.1)
## Analysis of Variance Table
##
## Model 1: raoQ ~ 1
## Model 2: raoQ ~ RDA1 + I(RDA1^2) + RDA2 + I(RDA2^2)
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 430 3181.0
## 2 426 2278.4 4 902.68 42.195 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(lm.rda.1)
##
## Call:
## lm(formula = raoQ ~ RDA1 + I(RDA1^2) + RDA2 + I(RDA2^2), data = arctic.cwm.rda.sites)
##
## Residuals:
## Min 1Q Median 3Q Max
## -9.6797 -1.4948 0.2389 1.4690 8.2247
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.979 0.147 142.681 < 2e-16 ***
## RDA1 -6.880 1.509 -4.559 6.72e-06 ***
## I(RDA1^2) -112.924 13.209 -8.549 2.25e-16 ***
## RDA2 -7.944 1.218 -6.523 1.95e-10 ***
## I(RDA2^2) -49.728 7.881 -6.310 7.00e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.313 on 426 degrees of freedom
## Multiple R-squared: 0.2838, Adjusted R-squared: 0.277
## F-statistic: 42.2 on 4 and 426 DF, p-value: < 2.2e-16
lm.cca.step<-stepAIC(lm.rda.1)
## Start: AIC=727.66
## raoQ ~ RDA1 + I(RDA1^2) + RDA2 + I(RDA2^2)
##
## Df Sum of Sq RSS AIC
## <none> 2278.4 727.66
## - RDA1 1 111.17 2389.5 746.19
## - I(RDA2^2) 1 212.95 2491.3 764.17
## - RDA2 1 227.59 2505.9 766.70
## - I(RDA1^2) 1 390.85 2669.2 793.90
raoQ.rda.1<-
ggplot()+
geom_smooth(data=arctic.cwm.rda.sites,
aes(x=RDA1, y=raoQ), formula = y~x,
method = "loess")+
geom_point(data=arctic.cwm.rda.sites,
aes(x=RDA1, y=raoQ, color=cluster))+
scale_color_manual(values=clus.cols.fuzzy)+
theme_bw()+
labs(x="RDA1", y="Rao's Quadratic Entropy")
raoQ.rda.1
raoQ.rda.2<-
ggplot()+
geom_smooth(data=arctic.cwm.rda.sites,
aes(x=RDA2, y=raoQ), formula = y~x,
method = "loess")+
geom_point(data=arctic.cwm.rda.sites,
aes(x=RDA2, y=raoQ, color=cluster))+
scale_color_manual(values=clus.cols.fuzzy)+
theme_bw()+
labs(x="RDA2", y="Rao's Quadratic Entropy")
raoQ.rda.2
Diversity metrics compared to environmental gradients for taxonomic
diversity
##taxonomy
arctic.cwm.cca.sites.tax$simpson<-as.numeric(arctic_env_fd$simpson)
lm.cca.0<-lm(simpson~ 1,
data=arctic.cwm.cca.sites.tax)
lm.cca.1<-lm(simpson~ CCA1 + I(CCA1^2) + CCA2 + I(CCA2^2),
data=arctic.cwm.cca.sites.tax)
anova(lm.cca.0, lm.cca.1)
## Analysis of Variance Table
##
## Model 1: simpson ~ 1
## Model 2: simpson ~ CCA1 + I(CCA1^2) + CCA2 + I(CCA2^2)
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 430 1.9641
## 2 426 1.5530 4 0.41117 28.197 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(lm.cca.1)
##
## Call:
## lm(formula = simpson ~ CCA1 + I(CCA1^2) + CCA2 + I(CCA2^2), data = arctic.cwm.cca.sites.tax)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.39112 -0.01927 0.01368 0.03797 0.11019
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.8759449 0.0040733 215.045 < 2e-16 ***
## CCA1 -0.0206280 0.0024860 -8.298 1.41e-15 ***
## I(CCA1^2) -0.0088179 0.0012681 -6.953 1.35e-11 ***
## CCA2 -0.0111958 0.0023405 -4.784 2.38e-06 ***
## I(CCA2^2) -0.0042692 0.0007753 -5.507 6.34e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.06038 on 426 degrees of freedom
## Multiple R-squared: 0.2093, Adjusted R-squared: 0.2019
## F-statistic: 28.2 on 4 and 426 DF, p-value: < 2.2e-16
lm.cca.step<-stepAIC(lm.cca.1)
## Start: AIC=-2414.78
## simpson ~ CCA1 + I(CCA1^2) + CCA2 + I(CCA2^2)
##
## Df Sum of Sq RSS AIC
## <none> 1.5530 -2414.8
## - CCA2 1 0.083416 1.6364 -2394.2
## - I(CCA2^2) 1 0.110539 1.6635 -2387.1
## - I(CCA1^2) 1 0.176258 1.7292 -2370.4
## - CCA1 1 0.251004 1.8040 -2352.2
summary(lm.cca.step)
##
## Call:
## lm(formula = simpson ~ CCA1 + I(CCA1^2) + CCA2 + I(CCA2^2), data = arctic.cwm.cca.sites.tax)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.39112 -0.01927 0.01368 0.03797 0.11019
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.8759449 0.0040733 215.045 < 2e-16 ***
## CCA1 -0.0206280 0.0024860 -8.298 1.41e-15 ***
## I(CCA1^2) -0.0088179 0.0012681 -6.953 1.35e-11 ***
## CCA2 -0.0111958 0.0023405 -4.784 2.38e-06 ***
## I(CCA2^2) -0.0042692 0.0007753 -5.507 6.34e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.06038 on 426 degrees of freedom
## Multiple R-squared: 0.2093, Adjusted R-squared: 0.2019
## F-statistic: 28.2 on 4 and 426 DF, p-value: < 2.2e-16
arctic.cwm.cca.sites.tax$CCA1<-as.numeric(arctic.cwm.cca.sites.tax$CCA1)
arctic.cwm.cca.sites.tax$CCA2<-as.numeric(arctic.cwm.cca.sites.tax$CCA2)
simpson.cca.1<-
ggplot()+
geom_smooth(data=arctic.cwm.cca.sites.tax,
aes(x=CCA1, y=simpson), formula = y~x,
method = "loess" )+
geom_point(data=arctic.cwm.cca.sites.tax,
aes(x=CCA1, y=simpson, color=cluster))+
scale_color_manual(values=clus.cols)+
theme_bw()+
labs(x="CCA1", y="Simpson Diversity Index")
simpson.cca.1
simpson.cca.2<-
ggplot()+
geom_smooth(data=arctic.cwm.cca.sites.tax,
aes(x=CCA2, y=simpson), formula = y~x,
method = "loess")+
geom_point(data=arctic.cwm.cca.sites.tax,
aes(x=CCA2, y=simpson, color=cluster))+
scale_color_manual(values=clus.cols)+
theme_bw()+
labs(x="CCA2", y="Simpson Diversity Index")
simpson.cca.2
Supplementary Figures
arctic_env_fd_2<-arctic_env_fd
arctic_env_fd_2$BottomTemp<-as.numeric(arctic_env_fd_2$BottomTemp)
arctic_env_fd_2$Year<-as.factor(arctic_env_fd_2$Year)
arctic_env_fd_melt<-melt(arctic_env_fd_2,
value.name = "value")
## Using Year, cluster, tax.cluster as id variables
View(arctic_env_fd_melt)
supp_fig_1<-ggplot(data=arctic_env_fd_melt,
aes(x=Year, y=value))+
geom_boxplot(aes(group=Year))+
facet_wrap(~variable, scales = "free", ncol = 2)
supp_fig_1