Effects of N fertilization on RMF, belowground biomass from Cleland et al. (2019)
Evan Perkowski
2022-11-29
Experiment leverages Nutrient Network study system to investigate how soil nutrient fertilization modifies belowground root mass fraction and total root biomass. Study quantifies RMF and total root biomass in 29 grassland sites in full factorial NPK combinations. Following analysis isolates Control (no nutrient addition) and +N (no P or K addition) plots:
Load packages
library(metafor)## Loading required package: Matrix## Loading required package: metadat##
## Loading the 'metafor' package (version 3.8-1). For an
## introduction to the package please type: help(metafor)library(MAd)
library(tidyverse)## ── Attaching packages
## ───────────────────────────────────────
## tidyverse 1.3.2 ──## ✔ ggplot2 3.3.6 ✔ purrr 0.3.5
## ✔ tibble 3.1.8 ✔ dplyr 1.0.10
## ✔ tidyr 1.2.1 ✔ stringr 1.4.1
## ✔ readr 2.1.3 ✔ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ tidyr::expand() masks Matrix::expand()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tidyr::pack() masks Matrix::pack()
## ✖ tidyr::unpack() masks Matrix::unpack()Data exploration
df <- read.csv("/Users/eaperkowski/Desktop/Cleland_2019_RMF.csv")
head(df)## site_code site_name collection_year elevation latitude longitude
## 1 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## 2 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## 3 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## 4 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## 5 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## 6 bldr.us Boulder South Campus 2011 1633 39.97202 -105.2335
## Aridity Light plot block treatment_name N P Kµ Rootsgperm2
## 1 0.366 0.7305256 4 1 NK 1 0 1 160.24611
## 2 0.366 0.7305256 18 2 NK 1 0 1 101.68877
## 3 0.366 0.7305256 16 2 PK 0 1 1 30.27665
## 4 0.366 0.7305256 7 1 NPK 1 1 1 77.71804
## 5 0.366 0.7305256 6 1 K 0 0 1 135.10583
## 6 0.366 0.7305256 8 1 NP 1 1 0 143.59417
## rootmassfraction
## 1 0.5140276
## 2 0.4649931
## 3 0.1169992
## 4 0.3450791
## 5 0.6960421
## 6 0.2885769Dataset includes full-factorial NPK additions. Subsetting to include only control and N addition plots:
df.summary <- df %>%
filter(treatment_name == "N" | treatment_name == "Control") %>%
group_by(site_code, treatment_name) %>%
summarize(x_rmf = mean(rootmassfraction, na.rm = TRUE),
sd_rmf = sd(rootmassfraction, na.rm = TRUE),
n_rmf = length(rootmassfraction),
x_rootbio = mean(Rootsgperm2, na.rm = TRUE),
sd_rootbio = sd(Rootsgperm2, na.rm = TRUE),
n_rootbio = length(Rootsgperm2)) %>%
pivot_wider(names_from = treatment_name, values_from = x_rmf:n_rootbio) %>%
rename(x_c_rmf = x_rmf_Control,
x_t_rmf = x_rmf_N,
sd_c_rmf = sd_rmf_Control,
sd_t_rmf = sd_rmf_N,
rep_c_rmf = n_rmf_Control,
rep_t_rmf = n_rmf_N,
x_c_rootbio = x_rootbio_Control,
x_t_rootbio = x_rootbio_N,
sd_c_rootbio = sd_rootbio_Control,
sd_t_rootbio = sd_rootbio_N,
rep_c_rootbio = n_rootbio_Control,
rep_t_rootbio = n_rootbio_N)
head(df.summary)## # A tibble: 6 × 13
## # Groups: site_code [6]
## site_code x_c_rmf x_t_rmf sd_c_rmf sd_t_rmf rep_c_…¹ rep_t…² x_c_r…³ x_t_r…⁴
## <chr> <dbl> <dbl> <dbl> <dbl> <int> <int> <dbl> <dbl>
## 1 bldr.us 0.455 0.532 0.400 0.251 2 2 87.4 231.
## 2 bnch.us 0.900 0.871 0.00985 0.0510 3 3 1264. 1154.
## 3 bogong.au 0.372 0.458 0.0290 0.0407 3 3 339. 819.
## 4 burrawan.au 0.429 0.227 0.0127 0.0796 3 3 136. 25.6
## 5 cbgb.us 0.271 0.283 0.0940 0.142 6 6 107. 103.
## 6 cdcr.us 0.457 0.358 0.248 0.0422 3 3 139. 137.
## # … with 4 more variables: sd_c_rootbio <dbl>, sd_t_rootbio <dbl>,
## # rep_c_rootbio <int>, rep_t_rootbio <int>, and abbreviated variable names
## # ¹rep_c_rmf, ²rep_t_rmf, ³x_c_rootbio, ⁴x_t_rootbioAnalyses
Root mass fraction:
df.rmf <- df.summary[1:7] %>%
rename(x_c = x_c_rmf,
x_t = x_t_rmf,
sd_c = sd_c_rmf,
sd_t = sd_t_rmf,
rep_c = rep_c_rmf,
rep_t = rep_t_rmf) %>%
metafor::escalc(measure = "ROM",
m1i = x_t, sd1i = sd_t, n1i = rep_t,
m2i = x_c, sd2i = sd_c, n2i = rep_c,
data = .,
append = TRUE,
var.names = c("logr", "logr_var")) %>%
mutate(logr_se = sqrt(logr_var) / sqrt(rep_t),
trait = "rmf")
head(df.rmf)##
## site_code x_c x_t sd_c sd_t rep_c rep_t logr
## 1 bldr.us 0.4547136 0.5322916 0.399977795 0.25066465 2 2 0.1575
## 2 bnch.us 0.9001407 0.8707471 0.009847249 0.05104679 3 3 -0.0332
## 3 bogong.au 0.3721667 0.4577132 0.028958297 0.04071909 3 3 0.2069
## 4 burrawan.au 0.4286752 0.2268936 0.012717256 0.07958051 3 3 -0.6362
## 5 cbgb.us 0.2709877 0.2834639 0.094028859 0.14205181 6 6 0.0450
## 6 cdcr.us 0.4571777 0.3578468 0.248107461 0.04224032 3 3 -0.2450
## logr_var logr_se trait
## 1 0.4978 0.49887461 rmf
## 2 0.0012 0.01987870 rmf
## 3 0.0047 0.03939633 rmf
## 4 0.0413 0.11733055 rmf
## 5 0.0619 0.10158860 rmf
## 6 0.1028 0.18512756 rmfRMF response ratio per site
## RMF response ratio per site
ggplot(data = df.rmf, aes(x = logr, y = site_code)) +
geom_point(size = 4) +
geom_vline(xintercept = 0, color = "red", linetype = "dashed") +
geom_errorbarh(aes(xmin = logr - logr_se, xmax = logr + logr_se)) +
labs(x = "Response to N fertilization", y = "Site") +
theme_bw(base_size = 12)
Root biomass per m^2:
df.rootbio <- df.summary[, c(1, 8:13)] %>%
rename(x_c = x_c_rootbio,
x_t = x_t_rootbio,
sd_c = sd_c_rootbio,
sd_t = sd_t_rootbio,
rep_c = rep_c_rootbio,
rep_t = rep_t_rootbio) %>%
metafor::escalc(
measure = "ROM",
m1i = x_t, sd1i = sd_t, n1i = rep_t,
m2i = x_c, sd2i = sd_c, n2i = rep_c,
data = .,
append = TRUE,
var.names = c("logr", "logr_var")) %>%
mutate(logr_se = sqrt(logr_var) / sqrt(rep_t),
trait = "rootbio")
head(df.rootbio)##
## site_code x_c x_t sd_c sd_t rep_c rep_t logr
## 1 bldr.us 87.43674 231.05094 82.99895 229.440573 2 2 0.9717
## 2 bnch.us 1263.60104 1153.58009 290.96178 156.654140 3 3 -0.0911
## 3 bogong.au 338.63679 819.43861 80.99709 229.865754 3 3 0.8837
## 4 burrawan.au 136.34385 25.55078 58.66165 4.048553 3 3 -1.6745
## 5 cbgb.us 106.93748 102.50577 36.45383 39.800163 6 6 -0.0423
## 6 cdcr.us 138.99092 137.41708 74.08509 61.753000 3 3 -0.0114
## logr_var logr_se trait
## 1 0.9436 0.68687273 rootbio
## 2 0.0238 0.08910833 rootbio
## 3 0.0453 0.12288166 rootbio
## 4 0.0701 0.15283252 rootbio
## 5 0.0445 0.08611377 rootbio
## 6 0.1620 0.23239268 rootbioRoot biomass per m^2 response ratio per site
## Root biomass response ratio per site
ggplot(data = df.rootbio, aes(x = logr, y = site_code)) +
geom_point(size = 4) +
geom_vline(xintercept = 0, color = "red", linetype = "dashed") +
geom_errorbarh(aes(xmin = logr - logr_se, xmax = logr + logr_se)) +
labs(x = "Response to N fertilization", y = "Site") +
theme_bw(base_size = 12)
Root biomass per m^2 response ratio across all sites
## Merge rootbio.all with rmf.all
cleland.total <- df.rootbio %>% full_join(df.rmf) %>%
dplyr::select(site_code, trait, everything())## Joining, by = c("site_code", "x_c", "x_t", "sd_c", "sd_t", "rep_c", "rep_t",
## "logr", "logr_var", "logr_se", "trait")## RMF response r
ggplot(data = cleland.total, aes(x = logr, y = trait)) +
geom_jitter(aes(size = 1/logr_se, fill = trait), height = 0.1, shape = 21) +
geom_boxplot(alpha = 0.6, width = 0.2, aes(fill = trait)) +
geom_vline(xintercept = 0, color = "red", linetype = "dashed") +
scale_y_discrete(labels = c("Root mass fraction",
expression("Root biomass (g m"^"-2"*")"))) +
labs(x = "Response to N fertilization", y = "Trait",
size = expression("Error"^"-1")) +
guides(fill = "none") +
theme_bw(base_size = 12)
Session info
sessionInfo()## R version 4.1.0 (2021-05-18)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] forcats_0.5.1 stringr_1.4.1 dplyr_1.0.10 purrr_0.3.5
## [5] readr_2.1.3 tidyr_1.2.1 tibble_3.1.8 ggplot2_3.3.6
## [9] tidyverse_1.3.2 MAd_0.8-3 metafor_3.8-1 metadat_1.2-0
## [13] Matrix_1.3-3
##
## loaded via a namespace (and not attached):
## [1] lubridate_1.9.0 lattice_0.20-44 assertthat_0.2.1
## [4] digest_0.6.29 utf8_1.2.2 R6_2.5.1
## [7] cellranger_1.1.0 backports_1.4.1 reprex_2.0.1
## [10] evaluate_0.15 highr_0.9 httr_1.4.4
## [13] pillar_1.8.1 rlang_1.0.6 googlesheets4_1.0.0
## [16] readxl_1.4.0 rstudioapi_0.13 jquerylib_0.1.4
## [19] rmarkdown_2.14 mathjaxr_1.6-0 labeling_0.4.2
## [22] googledrive_2.0.0 munsell_0.5.0 broom_0.8.0
## [25] compiler_4.1.0 modelr_0.1.8 xfun_0.30
## [28] pkgconfig_2.0.3 htmltools_0.5.2 tidyselect_1.2.0
## [31] fansi_1.0.3 crayon_1.5.2 withr_2.5.0
## [34] tzdb_0.3.0 dbplyr_2.1.1 grid_4.1.0
## [37] nlme_3.1-152 jsonlite_1.8.3 gtable_0.3.0
## [40] lifecycle_1.0.3 DBI_1.1.3 magrittr_2.0.3
## [43] scales_1.2.0 cli_3.4.1 stringi_1.7.8
## [46] farver_2.1.1 fs_1.5.2 xml2_1.3.3
## [49] bslib_0.3.1 ellipsis_0.3.2 generics_0.1.3
## [52] vctrs_0.5.1 tools_4.1.0 glue_1.6.2
## [55] hms_1.1.2 fastmap_1.1.0 yaml_2.3.5
## [58] timechange_0.1.1 colorspace_2.0-3 gargle_1.2.0
## [61] rvest_1.0.2 knitr_1.39 haven_2.5.0
## [64] sass_0.4.1