这是我第一次尝试在 R 中绘制 RDA。我的 RDA 数据是两个矩阵;一种与环境条件有关,例如 pH 值、温度等,另一种与不同的物种(每个站点每个物种的相对频率)有关。这些数据是在不同的站点获取的。我的响应变量是物种,我的解释变量是环境条件。
> chalut_large_mat #Species found at the stations
Actiniaria sp. Chionoecetes opilio Lycodes sp. Meganyctiphanes norvegica Pandalus borealis
EM14 0.04347826 0.1304348 0.04347826 0.04347826 0.6521739
EM16 0.00000000 0.0629275 0.00000000 0.00000000 0.8440492
BIC3 0.00000000 0.0000000 0.00000000 0.00000000 0.8431373
EM17 0.00000000 0.0000000 0.00000000 0.00000000 0.0000000
EM12 0.00000000 0.4076087 0.00000000 0.00000000 0.0000000
Pseudopleuronectes americanus Aspidophoroides monopterygius Crevette sp. Clymenella torquata Eteone sp.
EM14 0.08695652 0.00000000 0.00000000 0.0000000 0.00000000
EM16 0.00000000 0.03146375 0.00000000 0.0000000 0.00000000
BIC3 0.00000000 0.00000000 0.05882353 0.0000000 0.00000000
EM17 0.00000000 0.00000000 0.00000000 0.2142857 0.07142857
EM12 0.04347826 0.00000000 0.00000000 0.0000000 0.00000000
Mya arenaria Nuculana sp. Gammarus sp. Echinoidea Zoarcidae Sclerocrangon boreasÿ Pandalus montagui
EM14 0.0000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.0000000
EM16 0.0000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.0000000
BIC3 0.0000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.0000000
EM17 0.5714286 0.07142857 0.07142857 0.00000000 0.00000000 0.00000000 0.0000000
EM12 0.0000000 0.00000000 0.00000000 0.06521739 0.05978261 0.03804348 0.1521739
Sebastes sp.
EM14 0.0000000
EM16 0.0000000
BIC3 0.0000000
EM17 0.0000000
EM12 0.1684783
> chimie_large_mat #environmental conditions measured at every station
prof groupe lat long temp sal pH MPS chla phaeo oxygen_abs DO
BIC3 293 1 48.61018 -68.95022 5.9218 34.4832 7.885 11.29020 0.08353328 0.010738708 0.7734 53.06538
EM12 115 3 48.81703 -68.75740 3.3890 33.4904 7.986 15.95486 0.06623374 0.012311763 1.3739 135.49446
EM14 327 3 48.70100 -68.65102 6.2420 34.7703 7.919 15.72575 0.04919374 0.013765077 0.7829 49.63710
EM16 72 4 48.58262 -68.53828 7.4710 27.6365 8.197 9.96400 0.25115478 0.009154732 2.0803 284.46747
EM17 22 4 48.54998 -68.56083 6.6148 28.1667 8.207 11.71675 0.59121391 0.004247401 2.1196 270.89014
DIC AT phosphate NO2.NO3 pCO2
BIC3 2541.000 2275.107 2.481480 25.163658 1023.3
EM12 2380.742 2222.264 1.809933 18.547943 758.3
EM14 2494.135 2255.000 2.523323 24.916785 924.3
EM16 2297.125 2139.078 1.113962 7.995123 461.3
EM17 2221.665 2067.097 1.129026 10.197057 432.8
> ####ACR####
> acr<-rda(chalut_large_mat,chimie_large_mat,scale=FALSE)
> summary(acr)
Call:
rda(X = chalut_large_mat, Y = chimie_large_mat, scale = FALSE)
Partitioning of variance:
Inertia Proportion
Total 0.3105 1
Constrained 0.3105 1
Unconstrained 0.0000 0
Eigenvalues, and their contribution to the variance
Importance of components:
RDA1 RDA2 RDA3 RDA4
Eigenvalue 0.2307 0.07607 0.002925 0.0008345
Proportion Explained 0.7429 0.24499 0.009419 0.0026876
Cumulative Proportion 0.7429 0.98789 0.997312 1.0000000
Accumulated constrained eigenvalues
Importance of components:
RDA1 RDA2 RDA3 RDA4
Eigenvalue 0.2307 0.07607 0.002925 0.0008345
Proportion Explained 0.7429 0.24499 0.009419 0.0026876
Cumulative Proportion 0.7429 0.98789 0.997312 1.0000000
Scaling 2 for species and site scores
* Species are scaled proportional to eigenvalues
* Sites are unscaled: weighted dispersion equal on all dimensions
* General scaling constant of scores: 1.055675
Species scores
RDA1 RDA2 RDA3 RDA4
Actiniaria sp. -0.0094665 0.003555 0.0350172 -0.0053410
Chionoecetes opilio 0.1196407 0.296348 0.0175344 0.0253551
Lycodes sp. -0.0094665 0.003555 0.0350172 -0.0053410
Meganyctiphanes norvegica -0.0094665 0.003555 0.0350172 -0.0053410
Pandalus borealis -0.8158654 -0.104382 -0.0029265 0.0111654
Pseudopleuronectes americanus -0.0004895 0.038002 0.0621369 -0.0110708
Aspidophoroides monopterygius -0.0124340 -0.001963 -0.0067394 0.0225114
Crevette sp. -0.0237428 -0.007329 -0.0242368 -0.0357637
Clymenella torquata 0.1269299 -0.129706 0.0005282 0.0052062
Eteone sp. 0.0423100 -0.043235 0.0001761 0.0017354
Mya arenaria 0.3384797 -0.345882 0.0014086 0.0138833
Nuculana sp. 0.0423100 -0.043235 0.0001761 0.0017354
Gammarus sp. 0.0423100 -0.043235 0.0001761 0.0017354
Echinoidea 0.0276653 0.046338 -0.0118461 -0.0005831
Zoarcidae 0.0253599 0.042476 -0.0108590 -0.0005345
Sclerocrangon boreasÿ 0.0161381 0.027030 -0.0069102 -0.0003402
Pandalus montagui 0.0645524 0.108121 -0.0276410 -0.0013607
Sebastes sp. 0.0714687 0.119705 -0.0306025 -0.0015064
Site scores (weighted sums of species scores)
RDA1 RDA2 RDA3 RDA4
EM14 -0.2426 0.09112 0.897572 -0.136903
EM16 -0.4404 -0.06954 -0.238709 0.797357
BIC3 -0.4498 -0.13884 -0.459181 -0.677566
EM17 0.6601 -0.67457 0.002747 0.027076
EM12 0.4728 0.79183 -0.202429 -0.009965
Site constraints (linear combinations of constraining variables)
RDA1 RDA2 RDA3 RDA4
EM14 -0.2426 0.09112 0.897572 -0.136903
EM16 -0.4404 -0.06954 -0.238709 0.797357
BIC3 -0.4498 -0.13884 -0.459181 -0.677566
EM17 0.6601 -0.67457 0.002747 0.027076
EM12 0.4728 0.79183 -0.202429 -0.009965
Biplot scores for constraining variables
RDA1 RDA2 RDA3 RDA4
prof -0.7333 -0.20114 0.2821 -0.5850
groupe 0.6258 -0.02513 -0.7714 0.1125
lat -0.7837 -0.24602 -0.3465 0.4530
long 0.6651 -0.06890 -0.7343 -0.1173
我试过biplot
了,但我收到了这个错误消息,我知道这不是我想要表示我的数据的方式。
>biplot(acr, scaling = 1)
Error in biplot.rda(acr, scaling = 1, main = "RDA – scaling 1") :
'biplot.rda' not suitable for models with constraints
对于情节,我希望我的站点和我的物种都作为点,并将我的所有环境条件作为向量。我知道我可以用 ggplot 做到这一点,但我不知道怎么做。我试过了,但我只把我的站作为点,把我的物种作为矢量。我不知道如何将环境条件作为向量。这是我创建情节的代码,但显然是错误的。谢谢您的帮助!
acpl_summ<-summary(acr,scaling=1)
acpl_sites<-acpl_summ$sites[,1:2]
head(acpl_sites)
sites<-as_tibble(acpl_sites)%>%
mutate(station=rownames(acpl_sites))
sites
##Vector##
acpl_sp<-acpl_summ$species[,1:2]
acpl_sp
especes<-as_tibble(acpl_sp)%>%
mutate(espece=rownames(acpl_sp))
especes
##Plot##
acp_graph_sites<-ggplot(sites)+
geom_text(aes(x=RDA1,y=RDA2,label=station))+
xlab('RDA 1 (xx%)')+
ylab('RDA 2(xx%)')+
geom_hline(yintercept=0,linetype='dotted')+
geom_vline(xintercept=0,linetype='dotted')+
theme_bw()+
theme(panel.grid.major=element_blank(),
panel.grid.minor=element_blank(),
legend.position='top',
legend.text=element_text(size=7))+
scale_colour_discrete(name=NULL)
acp_graph_sites
exp_sp<-4 #expansion factor
acp_graph_sites_sp<-acp_graph_sites+
geom_text(data=especes,aes(x=RDA1/exp_sp,
y=RDA2/exp_sp,
label=espece))+
geom_segment(data=especes,
aes(x=0,y=0,xend=RDA1/exp_sp,yend=RDA2/exp_sp))
acp_graph_sites_sp