Birds_Inland_Sands_analysis_by_plot
Ron Chen
2024-01-29
## Loading required package: readxl## Loading required package: lubridate##
## Attaching package: 'lubridate'## The following objects are masked from 'package:data.table':
##
## hour, isoweek, mday, minute, month, quarter, second, wday, week,
## yday, year## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union## Loading required package: solartime## - Winter data excluded.
## - The counts of all passing / wintering species (not breeding in the sampled unit) were set to zero. Abreed_and_non_breed contains the non-breeders as well.
## - Species less likely to be interacting with sampling location were EXCLUDED
## - observations with counts greater than or equal to 10 standard deviations AND greater than or equal to 60 individuals were set to 1,
## under the assumption that these are migrating or local wandering phenomena.
## No outliers found in unit Arid South.
## 2 outliers set equal to 1 in Herbaceous and Dwarf-Shrub Vegetation unit.
## unit point_name year
## 1: Herbaceous and Dwarf-Shrub Vegetation Gamla Near 3 2014
## 2: Herbaceous and Dwarf-Shrub Vegetation Yavneel Far 2 2018
## SciName count_under_250 Z_score
## 1: Passer domesticus 62 11.90316
## 2: Passer hispaniolensis 60 11.35733
## No outliers found in unit Inland Sands.
## No outliers found in unit Loess Covered Areas in the Northern Negev.
## No outliers found in unit Mediterranean-Desert Transition Zone.
## 5 outliers set equal to 1 in Mediterranean Maquis unit.
## unit point_name year SciName
## 1: Mediterranean Maquis Beit Oren Near 1 2019 Curruca curruca
## 2: Mediterranean Maquis Beit Oren Near 3 2019 Curruca curruca
## 3: Mediterranean Maquis Kerem Maharal Near 3 2021 Chloris chloris
## 4: Mediterranean Maquis Kerem Maharal Near 3 2021 Garrulus glandarius
## 5: Mediterranean Maquis Nir Etzion Far 11 2021 Columba livia
## count_under_250 Z_score
## 1: 99 11.94535
## 2: 142 17.16927
## 3: 80 19.19266
## 4: 150 20.50138
## 5: 100 14.85118
## 2 outliers set equal to 1 in Negev Highlands unit.
## unit point_name year SciName
## 1: Negev Highlands Yeruham Far Slope 5 2018 Carpospiza brachydactyla
## 2: Negev Highlands Bislach Near Slope 6 2020 Passer hispaniolensis
## count_under_250 Z_score
## 1: 70 14.34802
## 2: 100 10.61343
## No outliers found in unit Planted Conifer Forest.Inland Sands Bird analysis
Monitoring started in 2017 and thus there are three campaigns until today: 2017, 2019 and 2021. There are three sites: Beer Milka, Secher and Shunra East. Factors are: - proximity to agriculture (Far vs. Near): only in Beer Milka. Other two sites are only Far. - dune status (Shifting vs. Semi-Shifting): in all three sites. There are a total of 12 plots, distributed as follows:
| site | agriculture | dunes | plot_num |
|---|---|---|---|
| Beer Milka | Far | Semi-Shifting | 2 |
| Beer Milka | Far | Shifting | 2 |
| Beer Milka | Near | Semi-Shifting | 2 |
| Beer Milka | Near | Shifting | 1 |
| Secher | Far | Semi-Shifting | 2 |
| Secher | Far | Shifting | 1 |
| Shunra East | Far | Semi-Shifting | 1 |
| Shunra East | Far | Shifting | 1 |
Raw data Total abundance: 720 Number of observations: 253 Total richness: 42
Filtered data Total abundance: 230 Number of observations: 94 Total richness: 23
ATTENTION: after filtration of rare species, there remain only 3 species. There are only 3 species in this unit with abundance > 10.
Model gma, abundance and richness
Richness done with rare species. Abundance and mean abundance done without rare species. Full models include cosine and sine of the time difference from June 21st, as well as time difference from sunrise (in radians).
richness
Explore data.
## [1] "RICHNESS WITH RARE SPECIES"
| richness | year_ct | site | agriculture | dunes | td_sc | cos_td_rad | sin_td_rad | h_from_sunrise | cos_hsun | sin_hsun | monitors_name | wind | precipitation | temperature | clouds | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Min. : 0.000 | Min. :5 | Beer Milka :21 | Far :27 | Semi-Shifting:21 | Min. :-1.263276 | Min. :0.1586 | Min. :-0.9873 | Length:36 | Min. :0.5686 | Min. :-0.07846 | Eyal Shochat :36 | 0 : 2 | 0 :12 | 0 : 0 | 0 : 0 | |
| 1st Qu.: 1.000 | 1st Qu.:5 | Secher : 9 | Near: 9 | Shifting :15 | 1st Qu.:-0.624197 | 1st Qu.:0.4325 | 1st Qu.:-0.9009 | Class :difftime | 1st Qu.:0.7808 | 1st Qu.: 0.13182 | Adi Domer : 0 | 1 : 5 | 3 : 0 | 1 : 5 | 1 : 0 | |
| Median : 1.500 | Median :7 | Shunra East: 6 | NA | NA | Median :-0.004195 | Median :0.6649 | Median :-0.7470 | Mode :numeric | Median :0.9509 | Median : 0.30898 | Asaf Mayrose : 0 | 2 : 4 | NA’s:24 | 2 : 6 | 2 : 0 | |
| Mean : 2.611 | Mean :7 | NA | NA | NA | Mean :-0.255375 | Mean :0.5599 | Mean :-0.7934 | NA | Mean :0.8892 | Mean : 0.34988 | Eliraz Dvir : 0 | 3 : 0 | NA | 3 : 1 | 3 : 0 | |
| 3rd Qu.: 3.000 | 3rd Qu.:9 | NA | NA | NA | 3rd Qu.: 0.129344 | 3rd Qu.:0.7079 | 3rd Qu.:-0.7049 | NA | 3rd Qu.:0.9912 | 3rd Qu.: 0.62473 | Eliraz Dvir and Yoav Barak: 0 | NA’s:25 | NA | NA’s:24 | NA’s:36 | |
| Max. :10.000 | Max. :9 | NA | NA | NA | Max. : 0.491806 | Max. :0.8140 | Max. :-0.5808 | NA | Max. :1.0000 | Max. : 0.82264 | Eran Banker : 0 | NA | NA | NA | NA | |
| NA | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA | (Other) : 0 | NA | NA | NA | NA |
observations missing for the 4 weather variables - exclude. Include sampling time of day.
## Warning in cor(P.anal[, lapply(X = .SD, FUN = as.numeric), .SDcols =
## IVs[1:12]], : the standard deviation is zero| richness | year_ct | site | agriculture | dunes | td_sc | cos_td_rad | sin_td_rad | h_from_sunrise | cos_hsun | sin_hsun | monitors_name | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| richness | 1.0000000 | -0.0135615 | -0.2600975 | 0.8055107 | -0.1160430 | 0.0905630 | 0.1018775 | 0.0596899 | -0.4964173 | 0.4116354 | -0.5104529 | NA |
| year_ct | -0.0135615 | 1.0000000 | 0.0000000 | 0.0000000 | 0.0000000 | -0.5863811 | -0.5680991 | -0.6061382 | -0.0753367 | 0.0715993 | -0.0732007 | NA |
| site | -0.2600975 | 0.0000000 | 1.0000000 | -0.4436070 | 0.0185535 | -0.2807880 | -0.2784000 | -0.2797564 | -0.2993936 | 0.2696945 | -0.3013447 | NA |
| agriculture | 0.8055107 | 0.0000000 | -0.4436070 | 1.0000000 | -0.0975900 | 0.2021769 | 0.2063717 | 0.1864087 | -0.2188171 | 0.1297906 | -0.2373704 | NA |
| dunes | -0.1160430 | 0.0000000 | 0.0185535 | -0.0975900 | 1.0000000 | -0.0050162 | 0.0008961 | -0.0171750 | -0.0738874 | 0.0499672 | -0.0770349 | NA |
| td_sc | 0.0905630 | -0.5863811 | -0.2807880 | 0.2021769 | -0.0050162 | 1.0000000 | 0.9974061 | 0.9853524 | 0.0569815 | -0.0612351 | 0.0533434 | NA |
| cos_td_rad | 0.1018775 | -0.5680991 | -0.2784000 | 0.2063717 | 0.0008961 | 0.9974061 | 1.0000000 | 0.9707262 | 0.0393742 | -0.0446465 | 0.0356772 | NA |
| sin_td_rad | 0.0596899 | -0.6061382 | -0.2797564 | 0.1864087 | -0.0171750 | 0.9853524 | 0.9707262 | 1.0000000 | 0.0981197 | -0.1007206 | 0.0944879 | NA |
| h_from_sunrise | -0.4964173 | -0.0753367 | -0.2993936 | -0.2188171 | -0.0738874 | 0.0569815 | 0.0393742 | 0.0981197 | 1.0000000 | -0.9605259 | 0.9983397 | NA |
| cos_hsun | 0.4116354 | 0.0715993 | 0.2696945 | 0.1297906 | 0.0499672 | -0.0612351 | -0.0446465 | -0.1007206 | -0.9605259 | 1.0000000 | -0.9434585 | NA |
| sin_hsun | -0.5104529 | -0.0732007 | -0.3013447 | -0.2373704 | -0.0770349 | 0.0533434 | 0.0356772 | 0.0944879 | 0.9983397 | -0.9434585 | 1.0000000 | NA |
| monitors_name | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA |
Fit Poisson glm, check for existence of overdispersion
## [1] "Estimating overdispersion parameter phi: (Res. Dev.)/(n-p) where n=number of observations; p=number of parameters in the model."## [1] "od = 0.419907633334411"Overdispersion parameter is < 1. Poisson more appropriate. Compare Poisson and negative binomial.
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached## df AIC
## m0.po 12 123.6344
## m0.nb 13 125.6355## [1] "poisson"
## [1] "neg bin"
negative binomial did not converge because there is underdispersion and poisson is better.
## boundary (singular) fit: see help('isSingular')mixed model is singular. Go with glm.
##
## Call:
## glm(formula = richness ~ agriculture * year_ct + dunes * year_ct +
## cos_td_rad + sin_td_rad + cos_hsun + sin_hsun + site, family = poisson,
## data = P.anal)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2304 -0.3246 -0.0637 0.3457 1.2679
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.81115 6.38362 -0.284 0.7766
## agricultureNear 2.64937 1.54878 1.711 0.0872 .
## year_ct 0.09271 0.16802 0.552 0.5811
## dunesShifting 1.89590 1.14298 1.659 0.0972 .
## cos_td_rad 2.44438 3.49098 0.700 0.4838
## sin_td_rad -3.00751 4.65376 -0.646 0.5181
## cos_hsun -1.23561 3.21472 -0.384 0.7007
## sin_hsun -2.82173 1.80813 -1.561 0.1186
## siteSecher 0.07287 0.41106 0.177 0.8593
## siteShunra East -0.41166 0.53002 -0.777 0.4373
## agricultureNear:year_ct -0.23537 0.24263 -0.970 0.3320
## year_ct:dunesShifting -0.30370 0.16337 -1.859 0.0630 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 73.183 on 35 degrees of freedom
## Residual deviance: 10.078 on 24 degrees of freedom
## AIC: 123.63
##
## Number of Fisher Scoring iterations: 5perform stepwise model selection of poisson model.
## Start: AIC=123.63
## richness ~ agriculture * year_ct + dunes * year_ct + cos_td_rad +
## sin_td_rad + cos_hsun + sin_hsun + site
##
## Df Deviance AIC
## - site 2 10.956 120.51
## - cos_hsun 1 10.225 121.78
## - sin_td_rad 1 10.503 122.06
## - cos_td_rad 1 10.576 122.13
## - agriculture:year_ct 1 11.011 122.57
## <none> 10.078 123.63
## - sin_hsun 1 12.539 124.10
## - year_ct:dunes 1 13.658 125.21
##
## Step: AIC=120.51
## richness ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## cos_hsun + sin_hsun + agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - cos_hsun 1 10.956 118.51
## - cos_td_rad 1 11.094 118.65
## - sin_td_rad 1 11.108 118.66
## - agriculture:year_ct 1 11.262 118.82
## - sin_hsun 1 12.594 120.15
## <none> 10.956 120.51
## - year_ct:dunes 1 13.886 121.44
##
## Step: AIC=118.51
## richness ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## sin_hsun + agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - cos_td_rad 1 11.094 116.65
## - sin_td_rad 1 11.108 116.66
## - agriculture:year_ct 1 11.263 116.82
## <none> 10.956 118.51
## - year_ct:dunes 1 13.908 119.46
## - sin_hsun 1 19.975 125.53
##
## Step: AIC=116.65
## richness ~ agriculture + year_ct + dunes + sin_td_rad + sin_hsun +
## agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - sin_td_rad 1 11.108 114.67
## - agriculture:year_ct 1 11.278 114.83
## <none> 11.094 116.65
## - year_ct:dunes 1 13.983 117.54
## - sin_hsun 1 19.984 123.54
##
## Step: AIC=114.66
## richness ~ agriculture + year_ct + dunes + sin_hsun + agriculture:year_ct +
## year_ct:dunes
##
## Df Deviance AIC
## - agriculture:year_ct 1 11.407 112.96
## <none> 11.108 114.67
## - year_ct:dunes 1 14.248 115.81
## - sin_hsun 1 22.224 123.78
##
## Step: AIC=112.96
## richness ~ agriculture + year_ct + dunes + sin_hsun + year_ct:dunes
##
## Df Deviance AIC
## <none> 11.407 112.96
## - year_ct:dunes 1 14.254 113.81
## - sin_hsun 1 24.766 124.32
## - agriculture 1 43.932 143.49Keep removing terms because \(\Delta AIC<2\). drop year:dunes
## Single term deletions
##
## Model:
## richness ~ agriculture + year_ct + dunes + sin_hsun
## Df Deviance AIC
## <none> 14.254 113.81
## agriculture 1 46.727 144.28
## year_ct 1 16.062 113.62
## dunes 1 14.734 112.29
## sin_hsun 1 25.683 123.24drop dunes.
## Single term deletions
##
## Model:
## richness ~ agriculture + year_ct + sin_hsun
## Df Deviance AIC
## <none> 14.734 112.29
## agriculture 1 49.098 144.66
## year_ct 1 16.594 112.15
## sin_hsun 1 25.811 121.37drop year.
## Single term deletions
##
## Model:
## richness ~ agriculture + sin_hsun
## Df Deviance AIC
## <none> 16.594 112.15
## agriculture 1 49.368 142.93
## sin_hsun 1 25.827 119.38Final model includes agriculture and sampling time of day. Final model:
##
## Call:
## glm(formula = richness ~ agriculture + sin_hsun, family = poisson,
## data = P.anal)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.87233 -0.27359 0.01039 0.33758 1.13786
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.8018 0.2082 3.851 0.000118 ***
## agricultureNear 1.2513 0.2196 5.698 1.21e-08 ***
## sin_hsun -1.2498 0.4313 -2.898 0.003757 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 73.183 on 35 degrees of freedom
## Residual deviance: 16.594 on 33 degrees of freedom
## AIC: 112.15
##
## Number of Fisher Scoring iterations: 4
| Observations | 36 |
| Dependent variable | richness |
| Type | Generalized linear model |
| Family | poisson |
| Link | log |
| χ²(2) | 56.590 |
| Pseudo-R² (Cragg-Uhler) | 0.801 |
| Pseudo-R² (McFadden) | 0.348 |
| AIC | 112.150 |
| BIC | 116.901 |
| exp(Est.) | 2.5% | 97.5% | z val. | p | |
|---|---|---|---|---|---|
| (Intercept) | 2.230 | 1.483 | 3.353 | 3.851 | 0.000 |
| agricultureNear | 3.495 | 2.272 | 5.374 | 5.698 | 0.000 |
| sin_hsun | 0.287 | 0.123 | 0.667 | -2.898 | 0.004 |
| Standard errors: MLE |
## Loading required package: extrafont## Registering fonts with R
There is a statistically significant effect for proximity to agriculture.
Near plots have on average 3.5921036 more species than far plots, which is 249.4767748 percent higher.
No statistically significant effect for the other factors tested.
geometric mean of abundance
ATTENTION: this is based on 3 species. There are only 3 species in this unit with abundance > 10.
Explore data. Include time of day.
## [1] "GEOMETRIC MEAN ABUNDANCE WITHOUT RARE SPECIES"## Warning: Removed 3 rows containing non-finite values (`stat_boxplot()`).
## Removed 3 rows containing non-finite values (`stat_boxplot()`).
## Removed 3 rows containing non-finite values (`stat_boxplot()`).
## Warning: Removed 3 rows containing non-finite values (`stat_boxplot()`).
## Removed 3 rows containing non-finite values (`stat_boxplot()`).
## Removed 3 rows containing non-finite values (`stat_boxplot()`).
| gma | year_ct | site | agriculture | td_sc | cos_td_rad | sin_td_rad | h_from_sunrise | cos_hsun | sin_hsun | monitors_name | wind | precipitation | temperature | clouds | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Min. : 1.000 | Min. :5 | Beer Milka :21 | Far :27 | Min. :-1.263276 | Min. :0.1586 | Min. :-0.9873 | Length:36 | Min. :0.5686 | Min. :-0.07846 | Eyal Shochat :36 | 0 : 2 | 0 :12 | 0 : 0 | 0 : 0 | |
| 1st Qu.: 2.000 | 1st Qu.:5 | Secher : 9 | Near: 9 | 1st Qu.:-0.624197 | 1st Qu.:0.4325 | 1st Qu.:-0.9009 | Class :difftime | 1st Qu.:0.7808 | 1st Qu.: 0.13182 | Adi Domer : 0 | 1 : 5 | 3 : 0 | 1 : 5 | 1 : 0 | |
| Median : 3.000 | Median :7 | Shunra East: 6 | NA | Median :-0.004195 | Median :0.6649 | Median :-0.7470 | Mode :numeric | Median :0.9509 | Median : 0.30898 | Asaf Mayrose : 0 | 2 : 4 | NA’s:24 | 2 : 6 | 2 : 0 | |
| Mean : 3.314 | Mean :7 | NA | NA | Mean :-0.255375 | Mean :0.5599 | Mean :-0.7934 | NA | Mean :0.8892 | Mean : 0.34988 | Eliraz Dvir : 0 | 3 : 0 | NA | 3 : 1 | 3 : 0 | |
| 3rd Qu.: 4.000 | 3rd Qu.:9 | NA | NA | 3rd Qu.: 0.129344 | 3rd Qu.:0.7079 | 3rd Qu.:-0.7049 | NA | 3rd Qu.:0.9912 | 3rd Qu.: 0.62473 | Eliraz Dvir and Yoav Barak: 0 | NA’s:25 | NA | NA’s:24 | NA’s:36 | |
| Max. :19.000 | Max. :9 | NA | NA | Max. : 0.491806 | Max. :0.8140 | Max. :-0.5808 | NA | Max. :1.0000 | Max. : 0.82264 | Eran Banker : 0 | NA | NA | NA | NA | |
| NA’s :3 | NA | NA | NA | NA | NA | NA | NA | NA | NA | (Other) : 0 | NA | NA | NA | NA |
extreme outlier where gma=19. Examine:
## unit subunit site year year_ct settlements agriculture
## 1: Inland Sands <NA> Shunra East 2017 5 <NA> Far
## habitat dunes land_use point_name date
## 1: <NA> Shifting <NA> Shunra East Far Shifting 1 2017-05-12
## datetime td_sc td_rad cos_td_rad sin_td_rad
## 1: 2017-05-12 06:15:00 0.3010365 -0.7057825 0.7611043 -0.6486296
## timediff_Jun21 monitors_name wind precipitation temperature clouds
## 1: -41 days Eyal Shochat <NA> <NA> <NA> <NA>
## h_from_sunrise cos_hsun sin_hsun pilot richness abundance gma
## 1: 0.36 hours 0.995562 0.09410831 FALSE 1 19 19Exclude this plot. Fit glm, compare gamma, gaussian (poisson inappropriate because response is not discrete)
Gamma seems better than gaussian. Fit fixed and mixed models.
## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
## Model failed to converge with max|grad| = 0.0530018 (tol = 0.002, component 1)Mixed model did not converge, use glm:
##
## Call:
## glm(formula = gma ~ agriculture * year_ct + dunes * year_ct +
## cos_td_rad + sin_td_rad + sin_hsun + cos_hsun + site, family = Gamma,
## data = P.anal)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.78874 -0.26688 -0.06482 0.15342 0.79350
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.38924 1.52848 -0.909 0.3742
## agricultureNear 0.39952 0.43647 0.915 0.3709
## year_ct 0.08400 0.04724 1.778 0.0906 .
## dunesShifting -0.12340 0.25074 -0.492 0.6280
## cos_td_rad 1.27164 0.81524 1.560 0.1345
## sin_td_rad -1.27885 1.06720 -1.198 0.2448
## sin_hsun -0.09136 0.48552 -0.188 0.8526
## cos_hsun -0.45039 0.96705 -0.466 0.6464
## siteSecher -0.16346 0.08088 -2.021 0.0569 .
## siteShunra East 0.10901 0.18358 0.594 0.5593
## agricultureNear:year_ct -0.08303 0.07004 -1.186 0.2497
## year_ct:dunesShifting 0.01756 0.03738 0.470 0.6436
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.2259429)
##
## Null deviance: 8.9535 on 31 degrees of freedom
## Residual deviance: 4.4842 on 20 degrees of freedom
## AIC: 112.15
##
## Number of Fisher Scoring iterations: 5perform stepwise model selection of Gamma model.
## Start: AIC=112.15
## gma ~ agriculture * year_ct + dunes * year_ct + cos_td_rad +
## sin_td_rad + sin_hsun + cos_hsun + site
##
## Df Deviance AIC
## - sin_hsun 1 4.4922 110.19
## - cos_hsun 1 4.5336 110.37
## - year_ct:dunes 1 4.5342 110.37
## - sin_td_rad 1 4.7998 111.55
## - agriculture:year_ct 1 4.8114 111.60
## <none> 4.4842 112.15
## - cos_td_rad 1 5.0302 112.57
## - site 2 5.8938 114.39
##
## Step: AIC=110.21
## gma ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## cos_hsun + site + agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - year_ct:dunes 1 4.5438 108.45
## - cos_hsun 1 4.6653 109.01
## - sin_td_rad 1 4.8067 109.66
## - agriculture:year_ct 1 4.8460 109.84
## <none> 4.4922 110.21
## - cos_td_rad 1 5.0731 110.89
## - site 2 6.2437 114.28
##
## Step: AIC=108.58
## gma ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## cos_hsun + site + agriculture:year_ct
##
## Df Deviance AIC
## - dunes 1 4.5459 106.59
## - cos_hsun 1 4.7207 107.42
## - sin_td_rad 1 4.8596 108.07
## - agriculture:year_ct 1 4.9281 108.39
## <none> 4.5438 108.58
## - cos_td_rad 1 5.1363 109.37
## - site 2 6.3032 112.87
##
## Step: AIC=106.6
## gma ~ agriculture + year_ct + cos_td_rad + sin_td_rad + cos_hsun +
## site + agriculture:year_ct
##
## Df Deviance AIC
## - cos_hsun 1 4.7210 105.47
## - sin_td_rad 1 4.8609 106.16
## - agriculture:year_ct 1 4.9414 106.56
## <none> 4.5459 106.60
## - cos_td_rad 1 5.1404 107.54
## - site 2 6.3068 111.32
##
## Step: AIC=105.84
## gma ~ agriculture + year_ct + cos_td_rad + sin_td_rad + site +
## agriculture:year_ct
##
## Df Deviance AIC
## - sin_td_rad 1 5.0904 105.70
## <none> 4.7210 105.84
## - cos_td_rad 1 5.4430 107.49
## - agriculture:year_ct 1 5.5936 108.25
## - site 2 6.4019 110.33
##
## Step: AIC=106.31
## gma ~ agriculture + year_ct + cos_td_rad + site + agriculture:year_ctkeep removing terms, because \(\Delta AIC<2\).
## Single term deletions
##
## Model:
## gma ~ agriculture + year_ct + cos_td_rad + site + agriculture:year_ct
## Df Deviance AIC
## <none> 5.0904 106.31
## cos_td_rad 1 5.8130 108.01
## site 2 6.7589 110.85
## agriculture:year_ct 1 5.6229 107.03drop agriculture:year
## Single term deletions
##
## Model:
## gma ~ agriculture + year_ct + site + cos_td_rad
## Df Deviance AIC
## <none> 5.6229 107.58
## agriculture 1 6.2733 108.57
## year_ct 1 6.3521 108.93
## site 2 7.4780 112.10
## cos_td_rad 1 6.0061 107.34drop cosine.
## Single term deletions
##
## Model:
## gma ~ agriculture + year_ct + site
## Df Deviance AIC
## <none> 6.0061 107.76
## agriculture 1 6.5819 108.25
## year_ct 1 6.3521 107.25
## site 2 8.4372 114.27drop year.
## Single term deletions
##
## Model:
## gma ~ agriculture + site
## Df Deviance AIC
## <none> 6.3521 107.61
## agriculture 1 6.8996 107.99
## site 2 8.8825 114.62keep removing terms, because \(\Delta AIC<2\). drop agriculture.
## Single term deletions
##
## Model:
## gma ~ site
## Df Deviance AIC
## <none> 6.8996 108.34
## site 2 8.9535 114.19In the final model only site remains.
No significant temporal or spatial trends detected in GMA.
abundance
ATTENTION: this is based on 3 species. There are only 3 species in this unit with abundance > 10.
Explore data.
## [1] "ABUNDANCE WITHOUT RARE SPECIES"
| abundance | year_ct | site | agriculture | td_sc | cos_td_rad | sin_td_rad | h_from_sunrise | cos_hsun | sin_hsun | monitors_name | wind | precipitation | temperature | clouds | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Min. : 0.000 | Min. :5 | Beer Milka :21 | Far :27 | Min. :-1.263276 | Min. :0.1586 | Min. :-0.9873 | Length:36 | Min. :0.5686 | Min. :-0.07846 | Eyal Shochat :36 | 0 : 2 | 0 :12 | 0 : 0 | 0 : 0 | |
| 1st Qu.: 1.000 | 1st Qu.:5 | Secher : 9 | Near: 9 | 1st Qu.:-0.624197 | 1st Qu.:0.4325 | 1st Qu.:-0.9009 | Class :difftime | 1st Qu.:0.7808 | 1st Qu.: 0.13182 | Adi Domer : 0 | 1 : 5 | 3 : 0 | 1 : 5 | 1 : 0 | |
| Median : 3.000 | Median :7 | Shunra East: 6 | NA | Median :-0.004195 | Median :0.6649 | Median :-0.7470 | Mode :numeric | Median :0.9509 | Median : 0.30898 | Asaf Mayrose : 0 | 2 : 4 | NA’s:24 | 2 : 6 | 2 : 0 | |
| Mean : 4.306 | Mean :7 | NA | NA | Mean :-0.255375 | Mean :0.5599 | Mean :-0.7934 | NA | Mean :0.8892 | Mean : 0.34988 | Eliraz Dvir : 0 | 3 : 0 | NA | 3 : 1 | 3 : 0 | |
| 3rd Qu.: 6.000 | 3rd Qu.:9 | NA | NA | 3rd Qu.: 0.129344 | 3rd Qu.:0.7079 | 3rd Qu.:-0.7049 | NA | 3rd Qu.:0.9912 | 3rd Qu.: 0.62473 | Eliraz Dvir and Yoav Barak: 0 | NA’s:25 | NA | NA’s:24 | NA’s:36 | |
| Max. :19.000 | Max. :9 | NA | NA | Max. : 0.491806 | Max. :0.8140 | Max. :-0.5808 | NA | Max. :1.0000 | Max. : 0.82264 | Eran Banker : 0 | NA | NA | NA | NA | |
| NA | NA | NA | NA | NA | NA | NA | NA | NA | NA | (Other) : 0 | NA | NA | NA | NA |
PHI>1, hence choose negative binomial. Fit fixed and mixed models. Choose mixed model if possible, otherwise choose a model with fixed-effects only.
## boundary (singular) fit: see help('isSingular')mixed model is singular. go with fixed effects glm.
##
## Call:
## glm.nb(formula = abundance ~ agriculture * year_ct + dunes *
## year_ct + cos_td_rad + sin_td_rad + sin_hsun + cos_hsun +
## site, data = P.anal, init.theta = 11.2265943, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.08784 -0.84611 -0.05132 0.37415 1.95844
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.34027 5.45301 0.979 0.3274
## agricultureNear -0.71201 1.40626 -0.506 0.6126
## year_ct -0.20899 0.16563 -1.262 0.2070
## dunesShifting 2.10540 0.93567 2.250 0.0244 *
## cos_td_rad -3.41135 2.89140 -1.180 0.2381
## sin_td_rad 4.38195 3.87517 1.131 0.2581
## sin_hsun -0.08653 1.63319 -0.053 0.9577
## cos_hsun 2.27396 2.98346 0.762 0.4459
## siteSecher 0.55097 0.32621 1.689 0.0912 .
## siteShunra East 0.31295 0.46400 0.674 0.5000
## agricultureNear:year_ct 0.31044 0.22669 1.369 0.1709
## year_ct:dunesShifting -0.26170 0.13733 -1.906 0.0567 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(11.2266) family taken to be 1)
##
## Null deviance: 93.230 on 35 degrees of freedom
## Residual deviance: 34.052 on 24 degrees of freedom
## AIC: 173.11
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 11.23
## Std. Err.: 8.63
##
## 2 x log-likelihood: -147.111
Perform stepwise model selection of mixed model.
## Start: AIC=171.11
## abundance ~ agriculture * year_ct + dunes * year_ct + cos_td_rad +
## sin_td_rad + sin_hsun + cos_hsun + site
##
## Df Deviance AIC
## - sin_hsun 1 34.054 169.11
## - cos_hsun 1 34.618 169.68
## - site 2 36.943 170.00
## - sin_td_rad 1 35.316 170.38
## - cos_td_rad 1 35.443 170.50
## - agriculture:year_ct 1 35.986 171.04
## <none> 34.052 171.11
## - year_ct:dunes 1 37.742 172.80
##
## Step: AIC=169.11
## abundance ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## cos_hsun + site + agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - site 2 36.934 168.03
## - sin_td_rad 1 35.309 168.40
## - cos_td_rad 1 35.503 168.60
## <none> 34.018 169.11
## - agriculture:year_ct 1 36.617 169.71
## - year_ct:dunes 1 37.813 170.91
## - cos_hsun 1 38.304 171.40
##
## Step: AIC=167.89
## abundance ~ agriculture + year_ct + dunes + cos_td_rad + sin_td_rad +
## cos_hsun + agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - sin_td_rad 1 35.281 166.95
## - cos_td_rad 1 35.826 167.49
## <none> 34.224 167.89
## - agriculture:year_ct 1 37.352 169.02
## - year_ct:dunes 1 37.851 169.52
## - cos_hsun 1 39.682 171.35
##
## Step: AIC=166.94
## abundance ~ agriculture + year_ct + dunes + cos_td_rad + cos_hsun +
## agriculture:year_ct + year_ct:dunes
##
## Df Deviance AIC
## - cos_td_rad 1 35.759 165.98
## <none> 34.720 166.94
## - agriculture:year_ct 1 37.009 167.23
## - year_ct:dunes 1 38.299 168.52
## - cos_hsun 1 39.867 170.09
##
## Step: AIC=165.96
## abundance ~ agriculture + year_ct + dunes + cos_hsun + agriculture:year_ct +
## year_ct:dunes
##
## Df Deviance AIC
## - agriculture:year_ct 1 38.122 165.39
## <none> 36.697 165.96
## - year_ct:dunes 1 40.832 168.10
## - cos_hsun 1 43.664 170.93
##
## Step: AIC=165.37
## abundance ~ agriculture + year_ct + dunes + cos_hsun + year_ct:dunes
##
## Df Deviance AIC
## <none> 37.218 165.37
## - year_ct:dunes 1 43.088 169.24
## - cos_hsun 1 50.896 177.05
## - agriculture 1 60.136 186.29dunes * year, agriculture and sampling time of day remain. The final model:
##
## Call:
## glm.nb(formula = abundance ~ agriculture + year_ct + dunes +
## cos_hsun + year_ct:dunes, data = P.anal, init.theta = 8.147282706,
## link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.1553 -0.8027 -0.2187 0.4935 2.0619
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.74659 1.04685 -1.668 0.095230 .
## agricultureNear 1.05307 0.22080 4.769 1.85e-06 ***
## year_ct -0.05424 0.08662 -0.626 0.531234
## dunesShifting 2.45559 0.93882 2.616 0.008906 **
## cos_hsun 3.34022 0.95175 3.510 0.000449 ***
## year_ct:dunesShifting -0.31888 0.13607 -2.343 0.019107 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(8.1473) family taken to be 1)
##
## Null deviance: 84.774 on 35 degrees of freedom
## Residual deviance: 37.218 on 30 degrees of freedom
## AIC: 167.37
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 8.15
## Std. Err.: 5.47
##
## 2 x log-likelihood: -153.373
Interpretation of abundance model:
## Error in glm.control(...) :
## unused argument (family = list("Negative Binomial(8.1473)", "log", function (mu)
## log(mu), function (eta)
## pmax(exp(eta), .Machine$double.eps), function (mu)
## mu + mu^2/.Theta, function (y, mu, wt)
## 2 * wt * (y * log(pmax(1, y)/mu) - (y + .Theta) * log((y + .Theta)/(mu + .Theta))), function (y, n, mu, wt, dev)
## {
## term <- (y + .Theta) * log(mu + .Theta) - y * log(mu) + lgamma(y + 1) - .Theta * log(.Theta) + lgamma(.Theta) - lgamma(.Theta + y)
## 2 * sum(term * wt)
## }, function (eta)
## pmax(exp(eta), .Machine$double.eps), expression({
## if (any(y < 0)) stop("negative values not allowed for the negative binomial family")
## n <- rep(1, nobs)
## mustart <- y + (y == 0)/6
## }), function (mu)
## all(mu > 0), function (eta)
## TRUE, function (object, nsim)
## {
## ftd <- fitted(object)
## rnegbin(nsim * length(ftd), ftd, .Theta)
## }))| Observations | 36 |
| Dependent variable | abundance |
| Type | Generalized linear model |
| Family | Negative Binomial(8.1473) |
| Link | log |
| χ²(NA) | NA |
| Pseudo-R² (Cragg-Uhler) | NA |
| Pseudo-R² (McFadden) | NA |
| AIC | 167.373 |
| BIC | 178.458 |
| exp(Est.) | 2.5% | 97.5% | z val. | p | |
|---|---|---|---|---|---|
| (Intercept) | 0.174 | 0.022 | 1.357 | -1.668 | 0.095 |
| agricultureNear | 2.866 | 1.860 | 4.419 | 4.769 | 0.000 |
| year_ct | 0.947 | 0.799 | 1.122 | -0.626 | 0.531 |
| dunesShifting | 11.653 | 1.851 | 73.379 | 2.616 | 0.009 |
| cos_hsun | 28.225 | 4.370 | 182.293 | 3.510 | 0.000 |
| year_ct:dunesShifting | 0.727 | 0.557 | 0.949 | -2.343 | 0.019 |
| Standard errors: MLE |
Test for difference of temporal trend from zero for each dune status.
## dunes year_ct.trend SE df z.ratio p.value
## Semi-Shifting -0.0542 0.0866 Inf -0.626 0.5312
## Shifting -0.3731 0.1070 Inf -3.486 0.0010
##
## Results are averaged over the levels of: agriculture
## P value adjustment: fdr method for 2 tests## dunes year_ct.trend SE df asymp.LCL asymp.UCL
## 1: Semi-Shifting -0.05423661 0.08662278 Inf -0.2240141 0.1155409
## 2: Shifting -0.37311230 0.10704437 Inf -0.5829154 -0.1633092## [1] 0.6885879## [1] 0.947208The trend of shifting sands is significantly different from zero, semi-shifting is not.
Test for pairwise differences in average abundance among levels of dunes, apply FDR correction for multiple comparisons.
## contrast year_ct estimate SE df z.ratio p.value
## (Semi-Shifting) - Shifting 7 -0.223 0.22 Inf -1.017 0.3093
##
## Results are averaged over the levels of: agriculture
## Results are given on the log (not the response) scale.There is a significant difference in average abundance among the two levels of distance from agriculture - near plots have on average 4.3394806 more individuals than far plots, which is 186.6426588%.
No effect of proximity to agriculture on the temporal trend of abundance.
No significant difference in average abundance among the two dune statuses - shifting (4.922029) and semi-shifting (3.9363821).
The temporal trend of shifting sands is significantly different from zero (P=0.001), semi-shifting is not. The rate is as follows: shifting sands decrease 31.14121% per year.
community analysis using package MVabund
Not performed, because after filtration of rare species we are left with very few species.