R包vegan的Mantel tests
Mantel tests是確定兩組距離測度矩陣（而非兩組變量矩陣）之間相關(guān)性的相關(guān)性測試方法揍瑟，用于判斷一個矩陣中的樣本距離與另一矩陣中的樣本距離是否相關(guān)铣鹏。Mantel tests零假設(shè)為響應(yīng)變量矩陣中對象之間的距離與解釋變量矩陣不存在相關(guān)肠牲，如果結(jié)果中p值顯著，則拒絕零假設(shè)指蚜，即存在相關(guān)性肛循，隨著一個矩陣中樣本之間距離的增加（或減少），另一矩陣中對應(yīng)樣本之間的距離也增加（或減少）柿究。
此外，Mantel方法還可用于檢驗假設(shè)或模型黄选。在這種模型測試方法中笛求，一個矩陣包含響應(yīng)數(shù)據(jù)，另一個矩陣代表了要測試的先驗?zāi)Ｐ停z驗的備擇假設(shè)）糕簿，如果找到了重要的Mantel統(tǒng)計信息，它們將為模型提供一些支持狡孔。
本篇重點介紹Mantel tests確定相關(guān)性的方法懂诗。
例如在群落分析中，為了探索群落物種組成是否與環(huán)境相關(guān)苗膝，經(jīng)常使用到Mantel tests殃恒。
群落分析中通常存在兩組變量矩陣，樣方-物種多度矩陣和樣方-環(huán)境測量矩陣辱揭。首先根據(jù)兩組變量矩陣計算樣方間的相異（距離）矩陣离唐，即分別獲得通過物種豐度計算的樣方距離（通常為Bray-curtis測度），以及通過某幾種環(huán)境參數(shù)計算的樣方距離（通常為歐幾里得距離）问窃。有時也會使用樣方間真實的地理距離直接作為某種距離測度亥鬓。
[圖片上傳失敗...(image-79b5f7-1609989983608)]
之后使用兩組距離測度矩陣執(zhí)行Mantel tests，例如確定樣方之間的群落組成差異是否與樣方之間的溫度差異或樣方之間的物理距離相關(guān)域庇，或者說“共變”嵌戈。這些測試可用于解決環(huán)境是針對微生物群落的“選擇”，還是存在強烈的距離衰減模式听皿，表明存在擴散限制熟呛。這些通常是生物地理學(xué)研究中的重要問題。
vegan包的Mantel tests方法
本篇同樣以群落分析為例尉姨，簡介R包vegan的Mantel tests庵朝。
假設(shè)存在如下數(shù)據(jù)集。第1列是樣方名稱又厉，第2-5列為各樣方中的環(huán)境參數(shù)（即鹽度九府、溫度等），第6-7列為各樣方的緯度和經(jīng)度覆致，第8列及之后為各樣方中的物種及其豐度昔逗。我們期望通過Mantel tests，查看對于該數(shù)據(jù)集篷朵，作用于物種變化的最主要因素勾怒，是由環(huán)境引起的“選擇”婆排，還是由地理因素的擴散限制所致。
加載R包笔链，如上所述段只，首先計算兩組樣方距離測度，然后執(zhí)行Mantel tests鉴扫。
library(vegan)
#讀取上述數(shù)據(jù)集
df <- read.csv('Your_OTU_table.csv', header= TRUE)
##計算距離
#根據(jù)物種豐度數(shù)據(jù)赞枕，計算樣方間的 Bray-curtis 距離
abund <- df[ ,8:ncol(df)]
dist.abund <- vegdist(abund, method = 'bray')
#根據(jù)環(huán)境測量指標，計算樣方間的歐幾里得距離
#這里只選擇了其中的溫度指標坪创，期望關(guān)注物種變化與溫度的相關(guān)性
temp <- df$Temperature
dist.temp <- dist(temp, method = 'euclidean')
#如果期望關(guān)注多種環(huán)境的協(xié)同作用炕婶，就選擇一個環(huán)境子集，計算樣方間的歐幾里得距離
#例如使用 4 種環(huán)境數(shù)據(jù)莱预，但此時需要執(zhí)行數(shù)據(jù)標準化柠掂，以消除量綱差異
env <- df[ ,2:5]
scale.env <- scale(env, center = TRUE, scale = TRUE)
dist.env <- dist(scale.env, method = 'euclidean')
#根據(jù)經(jīng)緯度，計算樣方間實際的地理距離
geo <- data.frame(df$Longitude, df$Latitude)
d.geo <- distm(geo, fun = distHaversine) #library(geosphere)
dist.geo <- as.dist(d.geo)
##執(zhí)行 Mantel tests依沮，詳情 ?mantel涯贞，以下為 3 個示例
#物種豐度和溫度的相關(guān)性，以 spearman 相關(guān)系數(shù)為例危喉，9999 次置換檢驗顯著性（Mantel tests 基于隨機置換的方法獲取 p 值）
abund_temp <- mantel(dist.abund, dist.temp, method = 'spearman', permutations = 9999, na.rm = TRUE)
abund_temp
#物種豐度和地理距離的相關(guān)性宋渔，以 spearman 相關(guān)系數(shù)為例，9999 次置換檢驗顯著性
abund_geo <- mantel(dist.abund, dist.geo, method = 'spearman', permutations = 9999, na.rm = TRUE)
abund_geo
#物種豐度和 4 種環(huán)境組合的相關(guān)性辜限，以 spearman 相關(guān)系數(shù)為例皇拣，9999 次置換檢驗顯著性
abund_env <- mantel(dist.abund, dist.env, method = 'spearman', permutations = 9999, na.rm = TRUE)
abund_env
基于物種豐度的距離矩陣與基于溫度指標的距離矩陣之間有很強的相關(guān)性（Mantel statistic R: 0.667，p value = 1e-04）薄嫡。換句話說审磁，隨著樣方在溫度方面的差異逐漸增大，它們在物種組成方面的差異也越來越大岂座。
#物種豐度和溫度的相關(guān)性
> abund_temp
Mantel statistic based on Spearman's rank correlation rho
Call:
mantel(xdis = dist.abund, ydis = dist.temp, method = "spearman", permutations = 9999, na.rm = TRUE)
Mantel statistic r: 0.677
Significance: 1e-04
Upper quantiles of permutations (null model):
90% 95% 97.5% 99%
0.148 0.198 0.246 0.290
Permutation: free
Number of permutations: 9999
基于物種豐度的距離矩陣與樣方間的地理距離沒有顯著關(guān)系（Mantel statistic R: 0.138态蒂，p value = 0.052）。因此可知费什，對于該測試數(shù)據(jù)集钾恢，不存在物種豐度的距離衰減效應(yīng)。
#物種豐度和地理距離的相關(guān)性
> abund_geo
Mantel statistic based on Spearman's rank correlation rho
Call:
mantel(xdis = dist.abund, ydis = dist.geo, method = "spearman", permutations = 9999, na.rm = TRUE)
Mantel statistic r: 0.1379
Significance: 0.0525
Upper quantiles of permutations (null model):
90% 95% 97.5% 99%
0.107 0.140 0.170 0.204
Permutation: free
Number of permutations: 9999
同時對于4種環(huán)境變量組合鸳址，累積的環(huán)境因素與群落物種組成高度相關(guān)（Mantel statistic r: 0.686, p value = 1e-04）瘩蚪。
#物種豐度和 4 種環(huán)境組合的相關(guān)性
> abund_env
Call:
mantel(xdis = dist.abund, ydis = dist.env, method = "spearman", permutations = 9999, na.rm = TRUE)
Mantel statistic r: 0.6858
Significance: 1e-04
Upper quantiles of permutations (null model):
90% 95% 97.5% 99%
0.151 0.201 0.244 0.292
Permutation: free
Number of permutations: 9999
綜上結(jié)論，對于該數(shù)據(jù)集稿黍，與地理距離相比疹瘦，群落物種組成與環(huán)境參數(shù)的相關(guān)性更強。因此在該系統(tǒng)中巡球，主要發(fā)生環(huán)境對群落作出的“選擇”言沐，地理因素的擴散限制相對微弱邓嘹。
作圖觀測相關(guān)性的示例
最后不妨作圖觀測變量間的關(guān)系，加深對這種相關(guān)性的理解险胰。
library(ggplot2)
#某物種與溫度的相關(guān)性汹押，橫軸溫度，縱軸物種豐度起便，顏色表示樣方的緯度
xx = ggplot(df, aes(x = Temperature, y = Pelagibacteraceae.OTU_307744)) +
geom_smooth(method = 'lm', alpha = 0.2, colour = 'black') +
geom_point(aes(colour = Latitude), size = 4) +
labs(y = 'Pelagibacteraceae (OTU 307744) (%)', x = 'Temperature (C)') +
theme( axis.text.x = element_text(face = 'bold',colour = 'black', size = 12),
axis.text.y = element_text(face = 'bold', size = 11, colour = 'black'),
axis.title= element_text(face = 'bold', size = 14, colour = 'black'),
panel.background = element_blank(),
panel.border = element_rect(fill = NA, colour = 'black'),
legend.title = element_text(size =12, face = 'bold', colour = 'black'),
legend.text = element_text(size = 10, face = 'bold', colour = 'black')) +
scale_colour_continuous(high = 'navy', low = 'salmon')
xx
#對于圖中的線性回歸
fit <- lm(df$Temperature~df$Pelagibacteraceae.OTU_307744)
summary(fit)
Call:
lm(formula = df$Temperature ~ df$Pelagibacteraceae.OTU_307744)
Residuals:
Min 1Q Median 3Q Max
-2.2053 -0.9336 -0.5215 0.5028 3.8232
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.4082 0.4476 0.912 0.372
df$Pelagibacteraceae.OTU_307744 1.3008 0.1280 10.165 1.45e-09 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.634 on 21 degrees of freedom
Multiple R-squared: 0.8311, Adjusted R-squared: 0.823
F-statistic: 103.3 on 1 and 21 DF, p-value: 1.454e-09
#分面圖展示多組變量的相關(guān)性棚贾，橫軸溫度，縱軸為多個物種的豐度榆综，顏色表示樣方的緯度
library(reshape2)
otus <- df[ ,1:11]
otus_melt <- melt(otus, id = c('Station', 'Salinity', 'Temperature', 'Oxygen', 'Nitrate', 'Latitude', 'Longitude'))
xx <- ggplot(otus_melt, aes(x = Temperature, y = value)) +
facet_wrap(.~variable, scales = 'free_y') +
geom_smooth(method = 'lm', alpha = 0.2, colour = 'black') +
geom_point(aes(colour = Latitude), size = 4) +
labs(y = 'Relative Abundance (%)', x = 'Temperature (C)') +
theme( axis.text.x = element_text(face = 'bold',colour = 'black', size = 12),
axis.text.y = element_text(face = 'bold', size = 10, colour = 'black'),
axis.title= element_text(face = 'bold', size = 14, colour = 'black'),
panel.background = element_blank(),
panel.border = element_rect(fill = NA, colour = 'black'),
legend.title = element_text(size =12, face = 'bold', colour = 'black'),
legend.text = element_text(size = 10, face = 'bold', colour = 'black'),
legend.position = 'top', strip.background = element_rect(fill = 'grey90', colour = 'black'),
strip.text = element_text(size = 9, face = 'bold')) +
scale_colour_continuous(high = 'navy', low = 'salmon')
xx
上述主要展示的變量間相關(guān)性的散點圖妙痹。
接下來是對于距離測度間的相關(guān)性。
#將上文獲得的距離測度鼻疮，轉(zhuǎn)化為數(shù)據(jù)框怯伊，一一對應(yīng)起來
aa <- as.vector(dist.abund)
tt <- as.vector(dist.temp)
gg <- as.vector(dist.geo)
mat <- data.frame(aa, tt, gg)
#基于物種豐度的距離與基于溫度指標的距離之間的相關(guān)性散點圖，上文已知二者顯著相關(guān)陋守；同時顏色表示樣方間地理距離
mm <- ggplot(mat, aes(y = aa, x = tt)) +
geom_point(size = 4, alpha = 0.75, colour = "black",shape = 21, aes(fill = gg/1000)) +
geom_smooth(method = "lm", colour = "black", alpha = 0.2) +
labs(x = "Difference in Temperature (C)", y = "Bray-Curtis Dissimilarity", fill = "Physical Separation (km)") +
theme( axis.text.x = element_text(face = "bold",colour = "black", size = 12),
axis.text.y = element_text(face = "bold", size = 11, colour = "black"),
axis.title= element_text(face = "bold", size = 14, colour = "black"),
panel.background = element_blank(),
panel.border = element_rect(fill = NA, colour = "black"),
legend.position = "top",
legend.text = element_text(size = 10, face = "bold"),
legend.title = element_text(size = 11, face = "bold")) +
scale_fill_continuous(high = "navy", low = "skyblue")
mm
#基于物種豐度的距離與樣方間地理距離之間的相關(guān)性散點圖，上文已知二者無相關(guān)性
mm <- ggplot(mat, aes(y = aa, x = gg/1000)) +
geom_point(size = 3, alpha = 0.5) +
labs(x = "Physical separation (km)", y = "Bray-Curtis Dissimilarity") +
theme( axis.text.x = element_text(face = "bold",colour = "black", size = 12),
axis.text.y = element_text(face = "bold", size = 11, colour = "black"),
axis.title= element_text(face = "bold", size = 14, colour = "black"),
panel.background = element_blank(),
panel.border = element_rect(fill = NA, colour = "black"))
mm
參考資料
Mantel Test in R：https://jkzorz.github.io/2019/07/08/mantel-test.html
Appendix 4: Graphical Description of MantelsTest and ANOSIM：https://www.mfe.govt.nz/publications/environmental-reporting/new-zealand-marine-environment-classification-overview/append-2
The Mantel test：https://mb3is.megx.net/gustame/hypothesis-tests/the-mantel-test
本文分享自微信公眾號 - 小白魚的生統(tǒng)筆記（gh_5f751e893315）利赋，作者：鯉小白
原文出處及轉(zhuǎn)載信息見文內(nèi)詳細說明水评，如有侵權(quán)，請聯(lián)系 yunjia_community@tencent.com 刪除媚送。
原始發(fā)表時間：2019-12-12
