R包vegan的Manteltests

R包vegan的Mantel tests
Mantel tests是確定兩組距離測度矩陣(而非兩組變數矩陣)之間相關性的相關性測試方法,用於判斷一個矩陣中的樣本距離與另一矩陣中的樣本距離是否相關。Mantel tests零假設為響應變數矩陣中物件之間的距離與解釋變數矩陣不存在相關,如果結果中p值顯著,則拒絕零假設,即存在相關性,隨著一個矩陣中樣本之間距離的增加(或減少),另一矩陣中對應樣本之間的距離也增加(或減少)。
此外,Mantel方法還可用於檢驗假設或模型。在這種模型測試方法中,一個矩陣包含響應資料,另一個矩陣代表了要測試的先驗模型(檢驗的備擇假設),如果找到了重要的Mantel統計資訊,它們將為模型提供一些支援。
本篇重點介紹Mantel tests確定相關性的方法。
例如在群落分析中,為了探索群落物種組成是否與環境相關,經常使用到Mantel tests
群落分析中通常存在兩組變數矩陣,樣方物種多度矩陣和樣方環境測量矩陣。首先根據兩組變數矩陣計算樣方間的相異(距離)矩陣,即分別獲得透過物種丰度計算的樣方距離(通常為Bray-curtis測度),以及透過某幾種環境引數計算的樣方距離(通常為歐幾里得距離)。有時也會使用樣方間真實的地理距離直接作為某種距離測度。

之後使用兩組距離測度矩陣執行Mantel tests,例如確定樣方之間的群落組成差異是否與樣方之間的溫度差異或樣方之間的物理距離相關,或者說“共變”。這些測試可用於解決環境是針對微生物群落的“選擇”,還是存在強烈的距離衰減模式,表明存在擴散限制。這些通常是生物地理學研究中的重要問題。
vegan包的Mantel tests方法
本篇同樣以群落分析為例,簡介RveganMantel tests
假設存在如下資料集。第1列是樣方名稱,第2-5列為各樣方中的環境引數(即鹽度、溫度等),第6-7列為各樣方的緯度和經度,第8列及之後為各樣方中的物種及其丰度。我們期望透過Mantel tests,檢視對於該資料集,作用於物種變化的最主要因素,是由環境引起的“選擇”,還是由地理因素的擴散限制所致。
載入R包,如上所述,首先計算兩組樣方距離測度,然後執行Mantel tests
library(vegan)
#讀取上述資料集

df <- read.csv('Your_OTU_table.csv', header= TRUE)
##計算距離

#根據物種丰度資料,計算樣方間的 Bray-curtis 距離

abund <- df[ ,8:ncol(df)]

dist.abund <- vegdist(abund, method = 'bray')
#根據環境測量指標,計算樣方間的歐幾里得距離

#這裡只選擇了其中的溫度指標,期望關注物種變化與溫度的相關性

temp <- df$Temperature

dist.temp <- dist(temp, method = 'euclidean')
#如果期望關注多種環境的協同作用,就選擇一個環境子集,計算樣方間的歐幾里得距離

#例如使用 4 種環境資料,但此時需要執行資料標準化,以消除量綱差異

env <- df[ ,2:5]

scale.env <- scale(env, center = TRUE, scale = TRUE)

dist.env <- dist(scale.env, method = 'euclidean')
#根據經緯度,計算樣方間實際的地理距離

geo <- data.frame(df$Longitude, df$Latitude)

d.geo <- distm(geo, fun = distHaversine)       #library(geosphere)

dist.geo <- as.dist(d.geo)
##執行 Mantel tests,詳情 ?mantel,以下為 3 個示例

#物種丰度和溫度的相關性,以 spearman 相關係數為例,9999 次置換檢驗顯著性(Mantel tests 基於隨機置換的方法獲取 p 值)

abund_temp <- mantel(dist.abund, dist.temp, method = 'spearman', permutations = 9999, na.rm = TRUE)

abund_temp
#物種丰度和地理距離的相關性,以 spearman 相關係數為例,9999 次置換檢驗顯著性

abund_geo <- mantel(dist.abund, dist.geo, method = 'spearman', permutations = 9999, na.rm = TRUE)

abund_geo
#物種丰度和 4 種環境組合的相關性,以 spearman 相關係數為例,9999 次置換檢驗顯著性

abund_env <- mantel(dist.abund, dist.env, method = 'spearman', permutations = 9999, na.rm = TRUE)

abund_env

基於物種丰度的距離矩陣與基於溫度指標的距離矩陣之間有很強的相關性(Mantel statistic R: 0.667p value = 1e-04)。換句話說,隨著樣方在溫度方面的差異逐漸增大,它們在物種組成方面的差異也越來越大。
#物種丰度和溫度的相關性

> 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

基於物種丰度的距離矩陣與樣方間的地理距離沒有顯著關係(Mantel statistic R: 0.138p value = 0.052)。因此可知,對於該測試資料集,不存在物種丰度的距離衰減效應。
#物種丰度和地理距離的相關性

> 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種環境變數組合,累積的環境因素與群落物種組成高度相關(Mantel statistic r: 0.686, p value = 1e-04)。
#物種丰度和 4 種環境組合的相關性

> 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

綜上結論,對於該資料集,與地理距離相比,群落物種組成與環境引數的相關性更強。因此在該系統中,主要發生環境對群落作出的“選擇”,地理因素的擴散限制相對微弱。
作圖觀測相關性的示例
最後不妨作圖觀測變數間的關係,加深對這種相關性的理解。
library(ggplot2)
#某物種與溫度的相關性,橫軸溫度,縱軸物種丰度,顏色表示樣方的緯度

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

#分面圖展示多組變數的相關性,橫軸溫度,縱軸為多個物種的丰度,顏色表示樣方的緯度

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

上述主要展示的變數間相關性的散點圖。
接下來是對於距離測度間的相關性。
#將上文獲得的距離測度,轉化為資料框,一一對應起來

aa <- as.vector(dist.abund)

tt <- as.vector(dist.temp)

gg <- as.vector(dist.geo)

mat <- data.frame(aa, tt, gg)
#基於物種丰度的距離與基於溫度指標的距離之間的相關性散點圖,上文已知二者顯著相關;同時顏色表示樣方間地理距離

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
#基於物種丰度的距離與樣方間地理距離之間的相關性散點圖,上文已知二者無相關性

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 Rhttps://jkzorz.github.io/2019/07/08/mantel-test.html
Appendix 4: Graphical Description of MantelsTest and ANOSIMhttps://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
連結

相關文章