之前介绍过拟合种面积关系(species–arearelationship, SAR)工具：
R——mmSAR对种面积关系进行拟合
今年3月份又出现了一个更强大的工具：sars
近期研究表明只使用单一的模型不能很好地拟合所有SAR数据，多个模型叠加可能更有实际意义。因此作者开发了sars：
可以使用线性或非线性的回归拟合20个不同的模型(mmSAR只有8个模型)；
还可以计算多个模型的平均曲线；
可用bootstrapping的方法得到置信区间。
SAR研究中使用最广泛的是幂律模型(power model)。但是一些研究已经发现大尺度上的SAR符合的是S型曲线(反曲型)。针对SAR模型不统一的情况，目前有两种策略，一是多个模型进行拟合，根据一定的标准选出效果最优(如AIC最小)的模型；二是多个模型拟合，取平均曲线。但是目前没有R包能实现。之前的两个包：
BAT可拟合三种SAR模型：线性、幂律和对数模型。
mmSAR可拟合8种模型，但是相比于目前超过20种的模型也不够用。
Sars相比于mmSAR的优势在于：
绘图更友好。可绘制加权的多模型曲线；
确定模型的形状(如线性、上凸、下凸、S型)；
是否有渐近线；
利用物种-地点丰度矩阵对Coleman’s (1981)的随机布局模型进行拟合；
建立了岛屿生物地理学的一般动态模型(general dynamic model, GDM)
20个模型
>install.packages("sars")
>library(sars)模型的拟合
非约束的Nelder–Mead最优化算法(R中的optim)最小化残差平方和来估计模型的参数。
代码语言：javascript复制##单个模型
>fit <- sar_loga(data = galap)
>summary(fit)
Model:
Logarithmic
Call:
S == c z * log(A)
Did the model converge: TRUE
Residuals:
0% 25% 50% 75% 100%
-164.000 -25.675 11.350 38.725 115.700
Parameters:
Estimate Std. Error t value Pr(>|t|) 2.5% 97.5%
c 0.2915341 34.5985767 0.0084262 0.9933958 -68.9056193 69.489
z 30.2802630 8.3203931 3.6392827 0.0026812 13.6394768 46.921
R-squared: 0.49, Adjusted R-squared: 0.41
AIC: 141.78, AICc: 143.78, BIC: 144.1
Observed shape: convex up, Asymptote: FALSE
Warning: The fitted values of the model contain negative values (i.e. negative species richness values)
>plot(fit)
#多个模型
>fitC <- sar_multi(data = galap, obj = c("power", "loga", "monod"))
Now attempting to fit the 3 SAR models:
-- multi_sars -------------------------------------------------- multi-model SAR --
> power : √
> loga : √
> monod : √
>plot(fitC)
模型的一般性和同质性检验>data(galap)
>fit <- sar_linear(data = galap, normaTest ="lillie", homoTest = "cor.fitted")
>summary(fit)
Model:
Linear model
Call:
S == c m*A
Did the model converge: TRUE
Residuals:
0% 25% 50% 75% 100%
-107.90 -52.10 -24.15 31.75 218.00
Parameters:
Estimate Std. Error t value Pr(>|t|) 2.5 % 97.5 %
c 70.314003 25.743130 2.731370 0.016228 15.100481 125.5275
m 0.185382 0.072884 2.543504 0.023410 0.029060 0.3417
R-squared: 0.32, Adjusted R-squared: 0.21
AIC: 146.36, AICc: 148.36, BIC: 148.68
Observed shape: linear, Asymptote: FALSE
Warning: The normality test selected indicated the model residuals are not normally distributed (i.e. P < 0.05)
>fit$normaTest
$test
[1] "lillie"
[[2]]
Lilliefors (Kolmogorov-Smirnov) normality test
data: res
D = 0.2146, p-value = 0.04725多模型曲线
将每一个成功拟合模型的预测丰度值乘以模型的权重(AIC,AICC,BIC等)，然后对所有模型的结果值求和，单个模型的线性组合构建多模型平均曲线。
代码语言：javascript复制>data(niering)
>fit <- sar_average(data= niering, obj =c("power","loga","koba","mmf","monod",
"negexpo","chapman","weibull3","asymp"),
normaTest = "none", homoTest = "none", neg_check = FALSE,
confInt = TRUE, ciN= 50)
Now attempting to fit the 9 SAR models:
-- multi_sars -------------------------------------------------- multi-model SAR --
> power : √
> loga : √
> koba : √
> mmf : √
> monod : √
> negexpo : √
> chapman : Warning: could not compute parameters statistics
> weibull3 : √
> asymp : x
1 models could not be fitted and have been excluded from the multi SAR
All models passed the model validation checks
8 remaining models used to construct the multi SAR:
Power, Logarithmic, Kobayashi, MMF, Monod, Negative exponential, Chapman Richards, Cumulative Weibull 3 par.
------------------------------------------------------------------------------------
Calculating sar_multi confidence intervals - this may take some time:
|==========================================================================| 100%
>par(mfrow = c(3,1))
#多个模型画图
>plot(fit, ModTitle = "a) Multimodel SAR")
#多模型平均曲线
>plot(fit, allCurves = FALSE, ModTitle =
"c) Multimodel SAR with confidence intervals", confInt = TRUE)
#每个模型所占权重
>plot(fit, type = "bar", ModTitle = "b) Model weights", cex.lab = 1.3)
