--- title: "HRRI: Reproduce All Six Paper Figures in R" author: "Mitra Ghotbi" output: rmarkdown::html_vignette: toc: true header-includes: | vignette: > %\VignetteIndexEntry{HRRI: Reproduce All Six Paper Figures in R} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include=FALSE} library(HRRI) knitr::opts_chunk$set(dev = "png", collapse=TRUE,comment="#>",fig.width=7.4, fig.height=4.8,dpi=200,out.width="100%",message=FALSE,warning=FALSE) library(HRRI) # Fast documentation build; export script defaults to all 2000 paper resamples. paper_n_boot <- 100L if (!exists("paper_n_boot", inherits = FALSE)) paper_n_boot <- 2000L library(HRRI) gal <- simulate_redox_holobiont(n_plot=2,n_depth=2,n_plant=3,n_time=40, seed=2026,scenario="flood_drain",disturbance_strength=.72,include_graph=TRUE, n_cycles=1L,disturbance_center=17, disturbance_width=5.5/sqrt(-2*log(.35))/40) fit <- function(sim) rri_pipeline_st(ROS_flux=sim$plant_data,Eh_stability=sim$Eh_stability, micro_data=log1p(sim$micro_gene_abundance),id=sim$id,time_col="time", group_cols=c("plot","depth","plant_id"),mode="snapshot", direction_anchor_phys="FvFm",direction_anchor_soil="Eh",direction_anchor_micro="mtrA") gal_res <- fit(gal) event_start <- min(gal$id$time[gal$id$phase=="disturbance"]) event_end <- max(gal$id$time[gal$id$phase=="disturbance"]) gal_scores <- attach_hrri_ids(gal_res$row_scores,gal$id) gal_scores$WFPS <- gal$forcing$WFPS gal_rec <- rri_recovery_metrics(gal_scores,time_col="time", group_cols=c("plot","depth","plant_id"),perturb_start=event_start, perturb_end=event_end,forcing_col="WFPS") props <- rri_property_scores(gal_res,rec=gal_rec,soil_df=gal$soil_data) acc_sim <- simulate_redox_holobiont(n_plot=24,n_depth=2,n_plant=3,n_time=40, seed=4096,scenario="flood_drain",disturbance_strength=.70,n_cycles=2L, disturbance_center=NULL,disturbance_width=.08) acc_res <- fit(acc_sim) acc <- rri_accuracy(acc_res$row_scores$RRI,acc_sim$latent_truth, cluster=acc_sim$id$plot,n_boot=paper_n_boot,n_perm=0,seed=20260913) paper_figures <- list( Fig1_framework=plot_rri_framework(), Fig2_identifiability=plot_rri_identifiability(), Fig3_trajectories=plot_rri_timeseries(gal,gal_res, perturb_start=event_start,perturb_end=event_end,forcing_threshold=.35,time_label="Time (days)"), Fig4_recovery_diagnostics=plot_rri_recovery_diagnostics(gal_res,gal$id,gal_rec, perturb_start=event_start,perturb_end=event_end), Fig5_operational_profile=plot_rri_properties(props,rec=gal_rec), Fig6_internal_agreement=plot_rri_accuracy(acc,base_size=9, score_label="Observation-derived score",target_label="Prescribed target", show_clusters=FALSE,cluster_label="Plots") ) show_figure <- function(p) { if(inherits(p,c("ggplot","patchwork"))) print(p) else for(q in p) print(q) } ``` # One figure-to-function map All six manuscript figures have native R implementations. Figures 1 and 2 are schematic and analytical illustrations; Figures 3–5 use the 12-trajectory gallery; Figure 6 uses a separate 24-plot simulation. These datasets must not be conflated. The scores and reference target share a synthetic generator: this is conditional internal agreement, not field validation. No Python installation is needed. | Paper figure | R function | |---|---| | 1: Complementary routes | `plot_rri_framework()` | | 2: Identifiability and accounting | `plot_rri_identifiability()` | | 3: Forcing and responses | `plot_rri_timeseries()` | | 4: Dynamics and availability | `plot_rri_recovery_diagnostics()` | | 5: Operational profile | `plot_rri_properties()` | | 6: Internal agreement | `plot_rri_accuracy()` | # Figure 1: complementary routes ```{r framework, fig.alt="Figure 1: complementary routes"} show_figure(paper_figures$Fig1_framework) ``` Declared inventories and exchange parameters define accessible electron capacity; reduced-species inventories and oxidation endpoints separately define oxygen demand. Observations yield domain scores and recovery diagnostics. The curve is schematic, and no inverse map from scores to reservoir parameters is implied. EAC and EDC remain separate directional capacities. # Figure 2: identifiability, sensitivity and accounting ```{r identifiability, fig.height=6.4, fig.alt="Figure 2: identifiability, sensitivity and accounting"} show_figure(paper_figures$Fig2_identifiability) ``` The one-pool illustration uses Q = 100, accessibility = 0.4, rate = 0.15 and duration = 10. One window leaves a curve of equivalent accessibility-rate pairs. Rate elasticity decreases near saturation; two-window inverse sensitivity quantifies local conditioning, not a confidence interval. The FeS comparison assumes oxidation of Fe(II) to Fe(III) and sulfide to sulfate. The equal-demand boundary is calculated from the declared stoichiometry, not fitted from scores. # Figure 3: forcing and responses ```{r trajectories, fig.height=6, fig.alt="Figure 3: forcing and responses"} show_figure(paper_figures$Fig3_trajectories) ``` The Gaussian forcing is shown explicitly with its 0.35 classification threshold. The highlighted days are `r event_start`–`r event_end`; nonzero forcing extends outside that classified interval. Eh, measured-inventory EAC and the score have separate axes and units. Missing observations remain gaps. EAC is not Cacc. # Figure 4: recovery diagnostics ```{r recovery, fig.width=8, fig.height=4.6, fig.alt="Figure 4: recovery diagnostics"} show_figure(paper_figures$Fig4_recovery_diagnostics) ``` Rows follow fixed identifier order and share a 0–1 score colour scale. Event markers are placed at the supplied observation times. Availability counts are computed from the recovery table. A finite value is not evidence of adequate precision, ecological validity or successful return. Rate and half-time share fit eligibility; `fit_status`, `n_fit` and the declared event window remain essential. # Figure 5: operational properties ```{r properties, fig.height=4.5, fig.alt="Figure 5: operational properties"} show_figure(paper_figures$Fig5_operational_profile) ``` The four descriptors are shown separately with provenance notes and numerical 0–1 ranges. These are not common physical units or a common favourable direction. The displayed timing and memory coverage counts come from the supplied recovery table, not manuscript constants. Inspect `props$property_table` for methods and `gal_rec` for trajectory-level support. No centre mean or polygon area is used. # Figure 6: conditional agreement ```{r agreement, fig.height=6.5, fig.alt="Figure 6: conditional agreement"} show_figure(paper_figures$Fig6_internal_agreement) ``` Small points are paired rows; open points are plot means. The score-on-target line and the identity line distinguish association from numerical agreement. Difference limits are the plot-mean bias plus or minus 1.96 SD of plot-mean differences; they are descriptive, not row-level limits or confidence intervals. Clustering follows independent plots, not individual plant-depth trajectories. Bootstrap draws condition on the fitted scores; the scoring pipeline is not refitted. This fast vignette uses **100 resamples**, so its bootstrap density and interval widths are illustrative. The supplied export script defaults to the paper's **2000 resamples**, with the same seed. The MSE partition uses population moments and closes algebraically. # Export without loss of quality Run explicitly in your R session (not executed during vignette builds): ```{r export, eval=FALSE} figure_dir <- file.path(getwd(), "HRRI_paper_figures") paper_n_boot <- 2000L source(system.file("examples", "export_paper_figures.R", package="HRRI"), local=TRUE) ``` PDF preserves vector geometry. Installing the suggested `svglite` package also provides editable-text SVG. PNG is exported at 600 dpi for previews; increasing raster resolution does not substitute for vector output. Use the supplied physical dimensions so text and panel spacing remain consistent. Scientific interpretation and provenance belong in the captions as well as the figures. # References ## Published methods Keiluweit, M., Wanzek, T., Kleber, M., Nico, P., & Fendorf, S. (2017). Anaerobic microsites have an unaccounted role in soil carbon stabilization. *Nature Communications*, **8**, 1771. Klüpfel, L., Piepenbrock, A., Kappler, A., & Sander, M. (2014). Humic substances as fully regenerable electron acceptors in recurrently anoxic environments. *Nature Geoscience*, **7**, 195–200. Kobayashi, K., & Salam, M. U. (2000). Comparing simulated and measured values using mean squared deviation and its components. *Agronomy Journal*, **92**, 345–352. Lin, L. I. (1989). A concordance correlation coefficient to evaluate reproducibility. *Biometrics*, **45**, 255–268. Sander, M., Hofstetter, T. B., & Gorski, C. A. (2015). Electrochemical analyses of redox-active iron minerals: a review of nonmediated and mediated approaches. *Environmental Science & Technology*, **49**, 5862–5878. Thompson, A., Chadwick, O. A., Rancourt, D. G., & Chorover, J. (2006). Iron-oxide crystallinity increases during soil redox oscillations. *Geochimica et Cosmochimica Acta*, **70**, 1710–1727. ## Companion manuscripts These describe the framework this package implements. None is published and two are under review; the entries are provisional and should be replaced with the published versions. Ghotbi, M., Ghotbi, M., Komluski, J., & Holtgrewe-Stukenbrock, E. H. HRRI: direction-aware diagnostics for soil–plant–microbiome redox recovery across hydroclimatic disturbances. *In preparation.* Ghotbi, M., Kolody, B. C., Ghotbi, M., & Holtgrewe-Stukenbrock, E. A Theory of Hydroclimatic Redox Resilience. *Submitted to Communications Earth & Environment.* — the source of the capacity, connectivity, kinetics and memory decomposition used throughout this vignette. Ghotbi, M., Ghotbi, M., Mühling, K. H., & Stukenbrock, E. H. Rhizosphere redox recovery after hydrological disturbances: mechanisms across the soil–plant–microbiome continuum. *Submitted to Soil Biology & Biochemistry.*