Method: Equivalent Soil Mass (ESM) – Haden, von
Haden, Yang & DeLucia (2020)
Compliance: VM0042 v2.2 / VCS v4.7
Both management zones processed: Project (WP) and
Control (BSL)
Key design features (v1.8):
Lower_cm values present in each
sheet — works for single-layer (0–30 cm) and multi-increment
data (e.g. 0–15, 15–30 cm or 0–10, 10–20, 20–30 cm) with no
configuration change required.UNC_fraction is applied uniformly when deducting from every
individual PAI’s area-weighted SOC stock. Per-PAI
UNC_fraction is still reported for transparency but is
not used in the deduction — the group value is.Output workbook (8 sheets):
| Sheet | VM0042 Reference | Contents |
|---|---|---|
README |
– | Version, method, assumptions, unit notes |
EQ_REF |
Eqs. 37-47, 74-79 | Equation reference for IQC/VVB auditing |
RUN_LOG |
– | Processing status per PAI x Year x Zone |
ALL_FD |
– | Combined fixed-depth data, both zones |
ALL_ESM |
– | Combined ESM-corrected data, both zones |
SUMMARY_SOC |
Eqs. 46/47 inputs | Mean / n / SD of ESM SOC (t CO2e/ha) per PAI x Year x Zone x Strata |
AREA_WT_SOC |
Eqs. 74, 46/47 | Area-weighted SOC per PAI x Year x Zone; PAI-level and group-level UNC both reported; group UNC used for deduction |
DELTA_NET_SOC |
Eqs. 37, 40, 44-47, 74 | Net SOC removal (CR) and reduction (ER) per PAI using group UNC deduction |
# ==============================================================================
# USER INPUTS -- edit only this chunk
# ==============================================================================
input_file_name <- "4842_SOC_data_for_ESM_V1_0.xlsx" ## << INPUT >>
# ESM reference depths are now AUTO-DETECTED from each sheet's Lower_cm values.
# Set esm_depths_override to NULL for auto-detection (recommended).
# To force fixed depths for ALL sheets set a numeric vector, e.g. c(0, 15, 30).
esm_depths_override <- NULL ## << INPUT >> NULL = auto-detect per sheet
min_core_length_cm <- 0 ## minimum core depth (cm)
extrapolation <- TRUE ## spline extrapolation?
year_col_name <- "Year" ## column holding sampling year
project_zone_value <- "project" ## scenario filter -- WP zone
control_zone_value <- "control" ## scenario filter -- BSL zone
# VM0042 Eq.74 parameters
# ASSUMPTION A1: t0.667 = 0.4307 (one-sided t at 66.7% CI, large-n limit, VM0042).
t0667 <- 0.4307 ## << INPUT >>
# VM0042 Eq.46/47: verification period denominator
# ASSUMPTION A2: verif_period_yrs derived from min(Year) and each monitoring year.
# Set to NA to auto-derive; set an integer (e.g. 5) to override.
verif_period_override <- NA ## << INPUT >> NA = auto, or integer override
# ==============================================================================
# FIXED CONSTANTS (VM0042 / IPCC / NIST)
# ==============================================================================
# ASSUMPTION A3: Strata_area values in the SOC file are in ACRES (US/CA practice).
ACRES_TO_HA <- 1 / 2.47105381 # NIST exact
# ASSUMPTION A4: ESM Cum_SOC_g_cm2 at max depth * 100 = t C/ha.
G_CM2_TO_T_HA <- 100
# ASSUMPTION A5: C to CO2e conversion = 44/12 (molar mass ratio CO2/C).
C_TO_CO2e <- 44 / 12
# ==============================================================================
# PATH RESOLUTION
# ==============================================================================
script_dir <- tryCatch(
dirname(rstudioapi::getSourceEditorContext()$path),
error = function(e) getwd()
)
input_file_path <- file.path(script_dir, input_file_name)
if (!file.exists(input_file_path))
stop("Cannot find '", input_file_name, "' in:\n ", script_dir)
run_stamp <- format(Sys.time(), "%Y%m%d_%H%M%S")
output_filename <- file.path(script_dir,
paste0("ESM_SOC_ALL_PAIs_", run_stamp, ".xlsx"))
cat("Input :", input_file_path, "\n")## Input : C:/Users/Quan/Carbon Friendly/Carbon Friendly - Verra/Project 4842 Sentry Ag Services Carbon Project/Validation/DATA/4842_SOC_data_for_ESM_V1_0.xlsx## Output: C:/Users/Quan/Carbon Friendly/Carbon Friendly - Verra/Project 4842 Sentry Ag Services Carbon Project/Validation/DATA/ESM_SOC_ALL_PAIs_20260402_185811.xlsx## Run : 2026-04-02 18:58:11 AESTcat("ESM depths: ", if (is.null(esm_depths_override)) "AUTO-DETECT per sheet" else
paste(esm_depths_override, collapse=", "), "\n")## ESM depths: AUTO-DETECT per sheetall_sheets <- readxl::excel_sheets(input_file_path)
pai_sheets <- all_sheets[grepl("PAI", all_sheets, ignore.case = TRUE)]
if (length(pai_sheets) == 0)
stop("No PAI sheets found in: ", input_file_name)
cat("PAI sheets:", paste(pai_sheets, collapse = ", "), "\n")## PAI sheets: PAI3, PAI5, PAI9, PAI11, PAI15, PAI16# ==============================================================================
# detect_esm_depths(): derive ESM reference depths from observed Lower_cm values
#
# The input file stores depths inverted: Lower_cm = surface (0), Upper_cm = bottom.
# After prepare_columns() swaps them the correct depths are in Lower_cm.
# This helper inspects the *post-swap* data frame and returns c(0, unique Lower_cm).
# Examples:
# Single-layer 0-30 cm -> c(0, 30)
# Two increments -> c(0, 15, 30)
# Three increments -> c(0, 10, 20, 30)
# ==============================================================================
detect_esm_depths <- function(df) {
# df has already been through prepare_columns() so Lower_cm is the TRUE bottom depth
bottoms <- sort(unique(na.omit(as.numeric(df[["Lower_cm"]]))))
bottoms <- bottoms[bottoms > 0]
if (length(bottoms) == 0) stop("No positive Lower_cm values found; cannot detect ESM depths.")
c(0, bottoms)
}
# ==============================================================================
# prepare_columns(): rename and depth-swap for V1_0 input file
#
# Column mapping (V1_0 file convention):
# Sample ID -> ID (sample identifier)
# Rep_ID -> Ref_ID (ESM reference profile)
# Rep_old/Rep.1 -> Rep (analytical replicate 1-6, only if present)
# Lower_cm (=0) -> Upper_cm (surface; file stores surface depth as Lower_cm=0)
# Upper_cm (=30) -> Lower_cm (core bottom; file stores bottom depth as Upper_cm=30)
#
# NOTE: 'Rep' is only removed when a renamed source (Rep_old / Rep.1) has already
# been promoted into that slot. The V1_0 file has a plain 'Rep' column which
# must be preserved as-is.
# ==============================================================================
prepare_columns <- function(df) {
nms <- names(df)
# Step 1: promote Rep_old or Rep.1 -> Rep (only if those columns exist)
rep_src <- if ("Rep_old" %in% nms) "Rep_old" else
if ("Rep.1" %in% nms) "Rep.1" else NA_character_
if (!is.na(rep_src)) {
nms[nms == rep_src] <- "Rep_ESM"
names(df) <- nms; nms <- names(df)
# Remove the original 'Rep' (FD rep) only when we have a replacement
if ("Rep" %in% nms) { df <- df[, nms != "Rep", drop = FALSE]; nms <- names(df) }
nms[nms == "Rep_ESM"] <- "Rep"
names(df) <- nms; nms <- names(df)
}
# If no Rep_old / Rep.1 exists, leave 'Rep' untouched
# Step 2: depth swap via temp name
# V1_0 stores: Lower_cm = 0 (surface), Upper_cm = 30 (bottom) -- inverted
# ESM expects: Upper_cm = 0 (surface), Lower_cm = 30 (bottom)
if (all(c("Upper_cm", "Lower_cm") %in% nms)) {
nms[nms == "Lower_cm"] <- "_uc_tmp"
nms[nms == "Upper_cm"] <- "Lower_cm"
nms[nms == "_uc_tmp"] <- "Upper_cm"
names(df) <- nms; nms <- names(df)
}
# Step 3: rename ID columns
nms[nms == "Sample ID"] <- "ID"
nms[nms == "Rep_ID"] <- "Ref_ID"
names(df) <- nms
df
}
# Safe NULL fallback: returns a if non-NULL and non-zero-length, else b.
# Deliberately avoids is.na() which returns logical(0) on lists/Style objects.
null_val <- function(a, b) if (!is.null(a) && length(a) > 0) a else b
# Spline: Hyman monotone spline with linear approx fallback
spline_safe <- function(x, y, xout) {
tryCatch(
spline(x = x, y = y, xout = xout, method = "hyman")$y,
error = function(e) approx(x = x, y = y, xout = xout, rule = 2)$y
)
}
# Pooled area-weighted variance of the mean:
# Var_aw_mean = SUM_s [ (A_s / A_tot)^2 * SD_s^2 / n_s ]
pooled_var_aw <- function(A, s2, n) {
A_tot <- sum(A, na.rm = TRUE)
if (A_tot == 0) return(NA_real_)
w <- A / A_tot
sum(w^2 * s2 / pmax(n, 1L), na.rm = TRUE)
}
# VM0042 Eq.74: UNC = min(1, sqrt(pooled_var) / |mean| * t0.667)
eq74_unc <- function(pooled_var, mean_soc, t = t0667) {
if (is.na(pooled_var) || is.na(mean_soc) || mean_soc == 0) return(NA_real_)
min(1, sqrt(pooled_var) / abs(mean_soc) * t)
}req_cols <- c("ID","Ref_ID","Rep","Upper_cm","Lower_cm",
"SOC_pct","SOM_pct","BD_g_cm3")
# Output containers
all_fd_rows <- list()
all_esm_rows <- list()
esm_samps <- list() # per-sample ESM SOC at max depth: both zones
strata_reg <- list() # pai_id -> data.frame(Strata, Strata_area_ha)
run_log <- list()
sheet_depth_log <- list() # record detected depths per sheet
for (sheet_name in pai_sheets) {
cat("\n\n---\n\n## PAI:", sheet_name, "\n\n")
# -- 1. Read -------------------------------------------------------------
raw_df <- tryCatch(
suppressMessages(readxl::read_excel(input_file_path, sheet = sheet_name)),
error = function(e) {
cat("::: {.warn-box}\nSKIPPED:", sheet_name, "--", conditionMessage(e), "\n:::\n")
NULL })
if (is.null(raw_df) || nrow(raw_df) == 0) {
cat("::: {.warn-box}\nSKIPPED: empty\n:::\n"); next }
names(raw_df) <- trimws(names(raw_df))
# -- 2. PAI ID -----------------------------------------------------------
pai_id <- if ("PAI ID" %in% names(raw_df)) {
v <- unique(na.omit(as.character(raw_df[["PAI ID"]])))
if (length(v) == 1L) v else sheet_name
} else sheet_name
cat("**PAI ID:**", pai_id, " \n")
# -- 3. Identify metadata columns BEFORE rename --------------------------
scen_col <- names(raw_df)[grepl("scenario|management.*zone",
tolower(names(raw_df)))][1]
strata_col <- names(raw_df)[grepl("^strata$|^strata_id$",
tolower(trimws(names(raw_df))))][1]
area_col <- names(raw_df)[grepl("strata_area", tolower(names(raw_df)))][1]
for (v in c("scen_col","strata_col","area_col")) {
val <- get(v)
if (length(val) == 0 || is.na(val)) assign(v, NA_character_)
}
cat("**Scenario col:**", scen_col,
"| **Strata col:**", strata_col,
"| **Area col:**", area_col, " \n")
# -- 4. Strata area registry (Project rows only, acres -> ha) ------------
if (!is.na(strata_col) && !is.na(area_col) && !is.na(scen_col)) {
src <- raw_df[grepl(paste0("^", project_zone_value, "$"),
tolower(trimws(as.character(raw_df[[scen_col]]))),
ignore.case = TRUE), ]
area_num <- suppressWarnings(as.numeric(src[[area_col]]))
sadf <- data.frame(Strata = as.character(src[[strata_col]]),
Strata_area_ha = area_num * ACRES_TO_HA,
stringsAsFactors = FALSE) |>
dplyr::group_by(Strata) |>
dplyr::summarise(Strata_area_ha = dplyr::first(Strata_area_ha), .groups = "drop") |>
dplyr::filter(!is.na(Strata_area_ha))
strata_reg[[pai_id]] <- sadf
cat("**Strata areas (ha):**",
paste(sadf$Strata, "=", round(sadf$Strata_area_ha, 1), collapse = " | "), " \n")
}
# -- 5. Strata lookup: Sample_ID -> Strata (pre-rename column names) -----
strata_lkp <- list()
if (!is.na(strata_col) && "Sample ID" %in% names(raw_df)) {
tmp <- raw_df[!duplicated(raw_df[["Sample ID"]]), c("Sample ID", strata_col)]
for (i in seq_len(nrow(tmp)))
strata_lkp[[as.character(tmp[["Sample ID"]][i])]] <- as.character(tmp[[strata_col]][i])
}
# -- 6. Prepare columns (rename + depth swap) ----------------------------
raw_df <- prepare_columns(raw_df)
# -- 7. Coerce numerics --------------------------------------------------
for (nc in c("Upper_cm","Lower_cm","SOC_pct","SOM_pct","BD_g_cm3","Rep")) {
if (nc %in% names(raw_df))
raw_df[[nc]] <- suppressWarnings(as.numeric(raw_df[[nc]]))
}
miss <- req_cols[!req_cols %in% names(raw_df)]
if (length(miss) > 0) {
cat("::: {.warn-box}\nSKIPPED -- missing cols:", paste(miss, collapse = ", "), "\n:::\n")
next
}
# -- 8. AUTO-DETECT ESM depths for this sheet ----------------------------
# Use override if supplied, otherwise detect from observed Lower_cm values.
# Detected depths = c(0, <all unique positive Lower_cm values>).
# Examples: single 0-30 layer -> c(0,30); two increments -> c(0,15,30), etc.
ESM_depths_cm <- if (!is.null(esm_depths_override)) {
esm_depths_override
} else {
detect_esm_depths(raw_df)
}
max_depth_cm <- max(ESM_depths_cm)
n_increments <- length(ESM_depths_cm) - 1L # number of depth intervals
cat("**ESM depths (cm):**", paste(ESM_depths_cm, collapse = " -> "),
paste0("(", n_increments, " increment", if(n_increments!=1) "s" else "", ")"), " \n")
sheet_depth_log[[sheet_name]] <- data.frame(
Sheet = sheet_name,
PAI_ID = pai_id,
ESM_depths = paste(ESM_depths_cm, collapse = ","),
N_increments = n_increments,
stringsAsFactors = FALSE
)
# Depth sanity check (Upper < Lower for all rows)
depth_ok <- all(raw_df$Upper_cm < raw_df$Lower_cm, na.rm = TRUE)
cat("**Depth check:** Upper =",
paste(sort(unique(na.omit(raw_df$Upper_cm))), collapse = ","),
"| Lower =",
paste(sort(unique(na.omit(raw_df$Lower_cm))), collapse = ","),
if (depth_ok) "-- OK \n" else "-- WARNING (Upper >= Lower found) \n")
# -- 9. Years ------------------------------------------------------------
if (year_col_name %in% names(raw_df)) {
raw_df[[year_col_name]] <- suppressWarnings(as.integer(raw_df[[year_col_name]]))
sampling_years <- sort(unique(na.omit(raw_df[[year_col_name]])))
} else {
sampling_years <- NA_integer_
}
cat("**Years:**", paste(sampling_years, collapse = ", "), " \n\n")
# -- 10. Year loop -------------------------------------------------------
for (yr in sampling_years) {
yr_label <- if (!is.na(yr)) as.character(yr) else "All"
year_df <- if (!is.na(yr))
raw_df[!is.na(raw_df[[year_col_name]]) & raw_df[[year_col_name]] == yr, ]
else raw_df
cat("### Year:", yr_label, "\n\n")
# -- 11. Zone loop (Project AND Control) --------------------------------
for (zone_val in c(project_zone_value, control_zone_value)) {
zone_label <- tools::toTitleCase(zone_val)
zone_df <- if (!is.na(scen_col) && scen_col %in% names(year_df)) {
year_df[grepl(paste0("^", zone_val, "$"),
tolower(trimws(as.character(year_df[[scen_col]]))),
ignore.case = TRUE), ]
} else {
if (zone_val == project_zone_value) year_df else data.frame()
}
if (nrow(zone_df) == 0) {
cat(" Zone:", zone_label, "-- no rows, skipped \n"); next
}
fd_core <- zone_df[, req_cols, drop = FALSE]
v_errs <- character(0)
if (anyNA(fd_core$ID) || anyNA(fd_core$Ref_ID) || anyNA(fd_core$Rep))
v_errs <- c(v_errs, "NA in ID/Ref_ID/Rep")
if (anyNA(fd_core$Upper_cm) || anyNA(fd_core$Lower_cm))
v_errs <- c(v_errs, "NA in depths")
if (length(v_errs) > 0) {
cat(" Zone:", zone_label, "-- skipped:", paste(v_errs, collapse = "; "), " \n")
run_log[[paste0(sheet_name,"_",yr_label,"_",zone_label)]] <- data.frame(
PAI_ID=pai_id, Sheet=sheet_name, Year=yr_label, Zone=zone_label,
Status="SKIPPED", Notes=paste(v_errs, collapse="; "),
Rows_in=nrow(zone_df), FD_records=0L, ESM_records=0L,
ESM_depths=paste(ESM_depths_cm, collapse=","),
stringsAsFactors=FALSE)
next
}
# -- Fixed-depth pre-processing -------------------------------------
fd_work <- fd_core
fd_work$Type <- "FD"
fd_work <- fd_work[, c("Type","ID","Rep","Ref_ID",
"Upper_cm","Lower_cm",
"SOC_pct","SOM_pct","BD_g_cm3")]
# Keep only profiles starting at surface (Upper_cm == 0)
filtered_FD <- fd_work |>
dplyr::group_by(ID, Rep) |>
dplyr::filter(min(Upper_cm, na.rm = TRUE) == 0) |>
dplyr::ungroup() |>
tidyr::drop_na() |>
dplyr::arrange(ID, Rep, Upper_cm, Lower_cm)
# Handle non-contiguous depth sequences: keep only the contiguous
# leading section up to the first gap
all_contig <- filtered_FD |>
dplyr::group_by(ID, Rep) |>
dplyr::filter(all(Upper_cm == dplyr::lag(Lower_cm, default = 0))) |>
dplyr::ungroup()
noncontig <- filtered_FD |>
dplyr::group_by(ID, Rep) |>
dplyr::filter(any(Upper_cm != dplyr::lag(Lower_cm, default = 0))) |>
dplyr::ungroup()
rmv_nc <- noncontig |>
dplyr::group_by(ID, Rep) |>
dplyr::filter(Upper_cm < Upper_cm[
which(Upper_cm != dplyr::lag(Lower_cm, default = 0))[1]]) |>
dplyr::ungroup()
filtered_FD <- dplyr::bind_rows(all_contig, rmv_nc) |>
dplyr::group_by(ID, Rep) |>
dplyr::filter(max(Lower_cm, na.rm = TRUE) >= min_core_length_cm) |>
dplyr::ungroup()
if (nrow(filtered_FD) == 0) {
cat(" Zone:", zone_label, "-- no rows after FD filtering \n"); next
}
# Prepend synthetic zero row for each (ID, Ref_ID, Rep)
zero_rows <- filtered_FD |>
dplyr::distinct(ID, Ref_ID, Rep) |>
dplyr::mutate(Type="FD", Upper_cm=0, Lower_cm=0,
SOC_pct=0, SOM_pct=0, BD_g_cm3=0)
modified_FD <- dplyr::bind_rows(zero_rows, filtered_FD) |>
dplyr::arrange(ID, Ref_ID, Rep, Upper_cm, Lower_cm) |>
dplyr::mutate(
Soil_g_cm2 = (Lower_cm - Upper_cm) * BD_g_cm3,
SOC_g_cm2 = (SOC_pct / 100) * Soil_g_cm2,
SOM_g_cm2 = (SOM_pct / 100) * Soil_g_cm2,
Min_Soil_g_cm2 = Soil_g_cm2 - SOM_g_cm2
)
cumulative_FD <- modified_FD |>
dplyr::group_by(ID, Rep) |>
dplyr::mutate(
Cum_Soil_g_cm2 = cumsum(Soil_g_cm2),
Cum_SOC_g_cm2 = cumsum(SOC_g_cm2),
Cum_SOM_g_cm2 = cumsum(SOM_g_cm2),
Cum_Min_Soil_g_cm2 = cumsum(Min_Soil_g_cm2)
) |>
dplyr::ungroup()
# -- ESM calculation ------------------------------------------------
# For each (sample, rep):
# 1. Build mean reference cumulative mass profile at ESM_depths_cm
# (average across all reps of the reference sample at each depth).
# 2. Spline-interpolate the sample's cumulative SOC / SOM onto those
# reference masses (mass-for-mass equivalence).
# 3. Store the ESM SOC at the deepest reference depth.
cumulative_ESM <- data.frame()
distinct_combos <- dplyr::distinct(cumulative_FD, ID, Ref_ID, Rep)
n_skip <- 0L
for (i in seq_len(nrow(distinct_combos))) {
cv <- distinct_combos[i, ]
crep <- dplyr::filter(cumulative_FD,
ID == cv$ID, Ref_ID == cv$Ref_ID, Rep == cv$Rep)
# Reference profile: mean cumulative min-soil mass across reps at each
# ESM depth node. Filter to nodes that are within the observed depth range.
cref_mean <- cumulative_FD |>
dplyr::filter(ID == cv$Ref_ID, Lower_cm %in% ESM_depths_cm) |>
dplyr::group_by(Upper_cm, Lower_cm) |>
dplyr::summarise(Cum_Min_Soil_g_cm2 = mean(Cum_Min_Soil_g_cm2, na.rm = TRUE),
.groups = "drop") |>
dplyr::mutate(Type="ESM", ID=cv$ID, Ref_ID=cv$Ref_ID, Rep=cv$Rep,
Soil_g_cm2=NA_real_, SOC_g_cm2=NA_real_,
SOM_g_cm2=NA_real_, Min_Soil_g_cm2=NA_real_,
BD_g_cm3=NA_real_, SOC_pct=NA_real_, SOM_pct=NA_real_,
Cum_Soil_g_cm2=NA_real_, Cum_SOC_g_cm2=NA_real_,
Cum_SOM_g_cm2=NA_real_)
# Restrict reference nodes to within the measured depth / mass range
cref_filt <- if (!extrapolation)
cref_mean[cref_mean$Cum_Min_Soil_g_cm2 <=
max(crep$Cum_Min_Soil_g_cm2, na.rm = TRUE), ]
else
cref_mean[cref_mean$Lower_cm <= max(crep$Lower_cm, na.rm = TRUE), ]
# Require at least 2 reference nodes and 2 sample points for interpolation
if (nrow(cref_filt) < 2 || nrow(crep) < 2) { n_skip <- n_skip + 1L; next }
# Spline: sample's cumulative mass (x) -> cumulative SOC/SOM (y)
# evaluated at each reference mass node
cref_filt$Cum_SOC_g_cm2 <- spline_safe(
crep$Cum_Min_Soil_g_cm2, crep$Cum_SOC_g_cm2, cref_filt$Cum_Min_Soil_g_cm2)
cref_filt$Cum_SOM_g_cm2 <- spline_safe(
crep$Cum_Min_Soil_g_cm2, crep$Cum_SOM_g_cm2, cref_filt$Cum_Min_Soil_g_cm2)
# Back-calculate layer properties from cumulative ESM values
cref_final <- cref_filt |>
dplyr::arrange(Lower_cm) |>
dplyr::mutate(
Min_Soil_g_cm2 = Cum_Min_Soil_g_cm2 - dplyr::lag(Cum_Min_Soil_g_cm2, default=0),
SOC_g_cm2 = Cum_SOC_g_cm2 - dplyr::lag(Cum_SOC_g_cm2, default=0),
SOM_g_cm2 = Cum_SOM_g_cm2 - dplyr::lag(Cum_SOM_g_cm2, default=0),
Soil_g_cm2 = Min_Soil_g_cm2 + SOM_g_cm2,
BD_g_cm3 = dplyr::if_else((Lower_cm - Upper_cm) > 0,
Soil_g_cm2 / (Lower_cm - Upper_cm), NA_real_),
SOC_pct = dplyr::if_else(Soil_g_cm2 > 0,
SOC_g_cm2 / Soil_g_cm2 * 100, NA_real_),
SOM_pct = dplyr::if_else(Soil_g_cm2 > 0,
SOM_g_cm2 / Soil_g_cm2 * 100, NA_real_),
Cum_Soil_g_cm2 = cumsum(Soil_g_cm2)
)
cumulative_ESM <- dplyr::bind_rows(cumulative_ESM, cref_final)
# Store per-sample ESM SOC at the deepest reference depth (max_depth_cm)
# ASSUMPTION A4: Cum_SOC_g_cm2 * G_CM2_TO_T_HA = t C/ha; * C_TO_CO2e = t CO2e/ha
esm_at_max <- cref_final[cref_final$Lower_cm == max_depth_cm, ]
if (nrow(esm_at_max) > 0) {
cum_soc <- esm_at_max$Cum_SOC_g_cm2[1]
esm_samps[[length(esm_samps) + 1L]] <- data.frame(
PAI_ID = pai_id,
Year = yr_label,
Zone = zone_label,
Sample_ID = as.character(cv$ID),
Ref_ID = as.character(cv$Ref_ID),
Rep = cv$Rep,
Strata = null_val(strata_lkp[[as.character(cv$ID)]], "Unknown"),
ESM_max_depth_cm = max_depth_cm,
Cum_SOC_g_cm2 = cum_soc,
SOC_tC_ha = cum_soc * G_CM2_TO_T_HA,
SOC_tCO2e_ha = cum_soc * G_CM2_TO_T_HA * C_TO_CO2e,
SOC_pct_ESM = esm_at_max$SOC_pct[1],
stringsAsFactors = FALSE
)
}
} # end ESM sample loop
# Post-processing: tag with PAI, Year, Zone; drop synthetic zero rows.
# FIX v1.8 (ALL_FD/ALL_ESM empty): openxlsx::writeData fails on dplyr
# tibbles carrying residual readxl attributes. Wrapping in as.data.frame()
# strips all dplyr/tibble class info so writeData writes rows correctly.
fd_out <- as.data.frame(
dplyr::mutate(
dplyr::filter(cumulative_FD, !(Upper_cm == 0 & Lower_cm == 0)),
PAI_ID = pai_id, Year = yr_label, Zone = zone_label, .before = 1
)
)
if (nrow(cumulative_ESM) > 0) {
esm_out <- as.data.frame(
dplyr::mutate(
dplyr::filter(cumulative_ESM, !(Upper_cm == 0 & Lower_cm == 0)),
PAI_ID = pai_id, Year = yr_label, Zone = zone_label, .before = 1
)
)
} else {
esm_out <- data.frame()
}
key <- paste0(sheet_name, "_", yr_label, "_", zone_label)
all_fd_rows[[key]] <- fd_out
all_esm_rows[[key]] <- esm_out
run_log[[key]] <- data.frame(
PAI_ID=pai_id, Sheet=sheet_name, Year=yr_label, Zone=zone_label,
Status="OK",
Notes=if (n_skip > 0) paste0(n_skip, " combos skipped (<2 ref/sample pts)") else "",
Rows_in=nrow(zone_df), FD_records=nrow(fd_out), ESM_records=nrow(esm_out),
ESM_depths=paste(ESM_depths_cm, collapse=","),
stringsAsFactors=FALSE)
cat(" Zone:", zone_label, "--",
nrow(fd_out), "FD |", nrow(esm_out), "ESM |",
"depths:", paste(ESM_depths_cm, collapse="-"),
if (n_skip > 0) paste0("| ", n_skip, " skipped") else "", " \n")
} # end zone loop
cat("\n")
} # end year loop
} # end PAI loopPAI ID: PAI3
Scenario col: Scenario | Strata col:
Strata | Area col: Strata_area
Strata areas (ha): Stratum 0 = 209.2 | Stratum 1 = 328
| Stratum 2 = 270.3
ESM depths (cm): 0 -> 1.03 -> 1.19 ->
1.31571428571429 -> 30 (4 increments)
Depth check: Upper = 0 | Lower =
1.03,1.19,1.31571428571429,30 – OK
Years: 2022, 2024
Zone: Project – 70 FD | 70 ESM | depths:
0-1.03-1.19-1.31571428571429-30
Zone: Control – 15 FD | 12 ESM | depths: 0-1.03-1.19-1.31571428571429-30
| 3 skipped
Zone: Project – 70 FD | 70 ESM | depths:
0-1.03-1.19-1.31571428571429-30
Zone: Control – 11 FD | 11 ESM | depths:
0-1.03-1.19-1.31571428571429-30
PAI ID: PAI5
Scenario col: Management zone | Strata
col: Strata_ID | Area col: Strata_area
Strata areas (ha): Stratum 0 = 87.6 | Stratum 1 = 56.3
| Stratum 2 = 105.3
ESM depths (cm): 0 -> 30 (1 increment)
Depth check: Upper = 0 | Lower = 30 – OK
Years: 2022, 2024
Zone: Project – 30 FD | 30 ESM | depths: 0-30
Zone: Control – 10 FD | 8 ESM | depths: 0-30 | 2 skipped
Zone: Project – 24 FD | 24 ESM | depths: 0-30
Zone: Control – 10 FD | 8 ESM | depths: 0-30 | 2 skipped
PAI ID: PAI9
Scenario col: Management zone | Strata
col: Strata_ID | Area col: Strata_area
Strata areas (ha): Stratum 0 = 79.8 | Stratum 1 = 97.6
| Stratum 2 = 14.1
ESM depths (cm): 0 -> 30 (1 increment)
Depth check: Upper = 0 | Lower = 30 – OK
Years: 2022, 2024
Zone: Project – 24 FD | 24 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
Zone: Project – 24 FD | 24 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
PAI ID: PAI11
Scenario col: Management zone | Strata
col: Strata_ID | Area col: Strata_area
Strata areas (ha): Stratum 0 = 132 | Stratum 1 = 100.2
| Stratum 2 = 148.9
ESM depths (cm): 0 -> 30 (1 increment)
Depth check: Upper = 0 | Lower = 30 – OK
Years: 2022, 2024
Zone: Project – 42 FD | 42 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
Zone: Project – 42 FD | 42 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
PAI ID: PAI15
Scenario col: Management zone | Strata
col: Strata_ID | Area col: Strata_area
Strata areas (ha): Stratum 0 = 83.3 | Stratum 1 = 858.4
| Stratum 2 = 1101.7
ESM depths (cm): 0 -> 30 (1 increment)
Depth check: Upper = 0 | Lower = 30 – OK
Years: 2022, 2024
Zone: Project – 210 FD | 210 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
Zone: Project – 275 FD | 275 ESM | depths: 0-30
Zone: Control – 6 FD | 4 ESM | depths: 0-30 | 2 skipped
PAI ID: PAI16
Scenario col: Management zone | Strata
col: Strata_ID | Area col: Strata_area
Strata areas (ha): Stratum 1 = 38.4 | Stratum 2 =
58.5
ESM depths (cm): 0 -> 30 (1 increment)
Depth check: Upper = 0 | Lower = 30 – OK
Years: 2022, 2024
Zone: Project – 12 FD | 12 ESM | depths: 0-30
Zone: Control – 8 FD | 4 ESM | depths: 0-30 | 4 skipped
Zone: Project – 12 FD | 12 ESM | depths: 0-30
Zone: Control – 8 FD | 4 ESM | depths: 0-30 | 4 skipped
# ==============================================================================
# Step 1: Combine all per-sample ESM records
# ==============================================================================
df_samps <- if (length(esm_samps) > 0) dplyr::bind_rows(esm_samps) else data.frame()
if (nrow(df_samps) == 0) stop("No ESM sample records produced. Check input data and RUN_LOG.")
# ==============================================================================
# Step 2: Average analytical reps per Sample_ID (anti-pseudoreplication)
# Each physical sample location contributes one value to strata statistics.
# ==============================================================================
df_by_samp <- df_samps |>
dplyr::group_by(PAI_ID, Year, Zone, Sample_ID, Strata, ESM_max_depth_cm) |>
dplyr::summarise(
n_reps = dplyr::n(),
SOC_tCO2e_ha = mean(SOC_tCO2e_ha, na.rm = TRUE),
SOC_tC_ha = mean(SOC_tC_ha, na.rm = TRUE),
SOC_pct_ESM = mean(SOC_pct_ESM, na.rm = TRUE),
.groups = "drop")
# ==============================================================================
# Step 3: Attach project-zone strata areas to sample-level data
# (Control samples receive NA for Strata_area_ha -- used in reporting only)
# ==============================================================================
if (length(strata_reg) > 0) {
area_lkp <- dplyr::bind_rows(lapply(names(strata_reg), function(p) {
d <- strata_reg[[p]]; d$PAI_ID <- p; d }))
df_by_samp <- dplyr::left_join(df_by_samp, area_lkp, by = c("PAI_ID","Strata"))
}
# ==============================================================================
# Step 4: SUMMARY_SOC -- mean / n / SD per PAI x Year x Zone x Strata
# ==============================================================================
summary_soc <- df_by_samp |>
dplyr::group_by(PAI_ID, Year, Zone, Strata, Strata_area_ha, ESM_max_depth_cm) |>
dplyr::summarise(
n_samples = dplyr::n(),
Mean_SOC_tCO2e_ha = mean(SOC_tCO2e_ha, na.rm = TRUE),
SD_SOC_tCO2e_ha = sd(SOC_tCO2e_ha, na.rm = TRUE),
SE_SOC_tCO2e_ha = SD_SOC_tCO2e_ha / sqrt(n_samples),
Mean_SOC_tC_ha = mean(SOC_tC_ha, na.rm = TRUE),
SD_SOC_tC_ha = sd(SOC_tC_ha, na.rm = TRUE),
Mean_SOC_pct = mean(SOC_pct_ESM, na.rm = TRUE),
SD_SOC_pct = sd(SOC_pct_ESM, na.rm = TRUE),
.groups = "drop"
) |>
dplyr::arrange(PAI_ID, Year, Zone, Strata)
cat("SUMMARY_SOC rows:", nrow(summary_soc), "\n")## SUMMARY_SOC rows: 62# ==============================================================================
# Step 5: AREA_WT_SOC -- area-weighted mean SOC per PAI x Year x Zone
#
# TWO uncertainty estimates are computed and reported:
#
# PAI_UNC_fraction : Eq.74 pooled within a single PAI's strata
# (for transparency / audit trail only)
#
# Group_UNC_fraction : Eq.74 pooled across ALL PAIs for the same Year x Zone
# (used for the actual deduction in DELTA_NET_SOC)
#
# ASSUMPTION A7: Using the group-level UNC means that uncertainty is assessed
# at the programme level rather than per-parcel, which is more conservative
# and consistent with pooled project-scale reporting under VM0042.
# ==============================================================================
# --- 5a. Per-PAI area-weighted SOC and PAI-level UNC -----------------------
pai_zone_wt <- summary_soc |>
dplyr::filter(!is.na(Strata_area_ha)) |>
dplyr::group_by(PAI_ID, Year, Zone) |>
dplyr::summarise(
ESM_max_depth_cm = dplyr::first(ESM_max_depth_cm),
Total_area_ha = sum(Strata_area_ha),
AreaWt_SOC_tCO2e_ha = sum(Mean_SOC_tCO2e_ha * Strata_area_ha) / sum(Strata_area_ha),
Pooled_Var_PAI = pooled_var_aw(
Strata_area_ha,
ifelse(is.na(SD_SOC_tCO2e_ha), 0, SD_SOC_tCO2e_ha)^2,
n_samples),
n_samples_total = sum(n_samples),
n_strata = dplyr::n(),
n_strata_no_SD = sum(is.na(SD_SOC_tCO2e_ha) | n_samples == 1L),
.groups = "drop"
) |>
dplyr::mutate(
PAI_UNC_fraction = mapply(eq74_unc, Pooled_Var_PAI, AreaWt_SOC_tCO2e_ha),
PAI_UNC_pct = PAI_UNC_fraction * 100,
Total_SOC_tCO2e = AreaWt_SOC_tCO2e_ha * Total_area_ha,
PAI_Notes = ifelse(n_strata_no_SD > 0,
paste0(n_strata_no_SD, " strata n=1 (SD=NA)"), "")
)
# --- 5b. Group (project-wide) UNC per Year x Zone --------------------------
# Pool variance across all PAIs treating each PAI as a "stratum" in the
# project-level area-weighted mean.
#
# Var_group = SUM_p [ (A_p / A_total)^2 * Var_PAI_p ]
# where Var_PAI_p is the within-PAI pooled variance of the area-weighted mean.
#
# Then: Group_UNC = min(1, sqrt(Var_group) / |mean_group_SOC| * t0.667)
group_unc_tbl <- pai_zone_wt |>
dplyr::group_by(Year, Zone) |>
dplyr::summarise(
Group_Total_area_ha = sum(Total_area_ha),
Group_AreaWt_SOC_tCO2e_ha = sum(AreaWt_SOC_tCO2e_ha * Total_area_ha) /
sum(Total_area_ha),
# Pool PAI-level variances weighted by area^2
Group_Pooled_Var = pooled_var_aw(
Total_area_ha,
Pooled_Var_PAI,
rep(1L, dplyr::n())),
n_PAIs = dplyr::n(),
.groups = "drop"
) |>
dplyr::mutate(
Group_UNC_fraction = mapply(eq74_unc, Group_Pooled_Var, Group_AreaWt_SOC_tCO2e_ha),
Group_UNC_pct = Group_UNC_fraction * 100
)
# --- 5c. Join group UNC back to per-PAI table --------------------------------
pai_zone_wt <- dplyr::left_join(
pai_zone_wt,
dplyr::select(group_unc_tbl, Year, Zone,
Group_UNC_fraction, Group_UNC_pct,
Group_Pooled_Var,
Group_Total_area_ha, Group_AreaWt_SOC_tCO2e_ha),
by = c("Year","Zone")
) |>
# Deduct using GROUP uncertainty (not PAI-level)
dplyr::mutate(
SOC_net_tCO2e_ha = AreaWt_SOC_tCO2e_ha *
(1 - ifelse(is.na(Group_UNC_fraction), 0, Group_UNC_fraction)),
SOC_net_tCO2e = SOC_net_tCO2e_ha * Total_area_ha
)
# --- 5d. PROJECT summary row (one per Year x Zone) --------------------------
proj_level <- pai_zone_wt |>
dplyr::group_by(Year, Zone) |>
dplyr::summarise(
PAI_ID = "PROJECT",
ESM_max_depth_cm = dplyr::first(ESM_max_depth_cm),
Total_area_ha = sum(Total_area_ha),
AreaWt_SOC_tCO2e_ha = dplyr::first(Group_AreaWt_SOC_tCO2e_ha),
Pooled_Var_PAI = dplyr::first(Group_Pooled_Var),
n_samples_total = sum(n_samples_total),
n_strata = sum(n_strata),
n_strata_no_SD = sum(n_strata_no_SD),
PAI_UNC_fraction = dplyr::first(Group_UNC_fraction),
PAI_UNC_pct = dplyr::first(Group_UNC_pct),
Total_SOC_tCO2e = AreaWt_SOC_tCO2e_ha * Total_area_ha,
PAI_Notes = "PROJECT row: group pooled",
Group_UNC_fraction = dplyr::first(Group_UNC_fraction),
Group_UNC_pct = dplyr::first(Group_UNC_pct),
Group_Total_area_ha = Total_area_ha,
Group_AreaWt_SOC_tCO2e_ha = AreaWt_SOC_tCO2e_ha,
SOC_net_tCO2e_ha = AreaWt_SOC_tCO2e_ha *
(1 - ifelse(is.na(dplyr::first(Group_UNC_fraction)), 0,
dplyr::first(Group_UNC_fraction))),
SOC_net_tCO2e = SOC_net_tCO2e_ha * Total_area_ha,
.groups = "drop"
)
area_wt_soc <- dplyr::bind_rows(pai_zone_wt, proj_level) |>
dplyr::arrange(Year, Zone, PAI_ID)
cat("AREA_WT_SOC rows:", nrow(area_wt_soc), "\n")## AREA_WT_SOC rows: 28# ==============================================================================
# Step 6: DELTA_NET_SOC -- WP minus BSL per PAI per Year
#
# Deduction uses Group_UNC_fraction (same value for all PAIs in same Year x Zone).
# ==============================================================================
delta_net_soc <- if (nrow(area_wt_soc) > 0) {
aw_pai <- area_wt_soc |>
dplyr::filter(PAI_ID != "PROJECT") |>
dplyr::select(PAI_ID, Year, Zone, Total_area_ha,
AreaWt_SOC_tCO2e_ha, Group_UNC_fraction, Total_SOC_tCO2e)
wp <- aw_pai |> dplyr::filter(Zone == "Project") |>
dplyr::rename(WP_area_ha = Total_area_ha,
WP_SOC_ha = AreaWt_SOC_tCO2e_ha,
WP_Group_UNC = Group_UNC_fraction,
WP_SOC_tCO2e = Total_SOC_tCO2e) |>
dplyr::select(-Zone)
bsl <- aw_pai |> dplyr::filter(Zone == "Control") |>
dplyr::rename(BSL_area_ha = Total_area_ha,
BSL_SOC_ha = AreaWt_SOC_tCO2e_ha,
BSL_Group_UNC = Group_UNC_fraction,
BSL_SOC_tCO2e = Total_SOC_tCO2e) |>
dplyr::select(-Zone)
delta <- dplyr::left_join(wp, bsl, by = c("PAI_ID","Year")) |>
dplyr::mutate(
Year_int = suppressWarnings(as.integer(Year)),
# -- Stock comparison at each sampled year (for reporting / audit) ------
Delta_SOC_ha = WP_SOC_ha - BSL_SOC_ha, # t CO2e/ha stock difference
Net_SOC_tCO2e = Delta_SOC_ha * WP_area_ha, # t CO2e total stock difference
# ASSUMPTION A6: Group_UNC_combined = max(WP_Group_UNC, BSL_Group_UNC)
# Taking the larger of the two zone UNCs is the conservative choice.
Group_UNC_combined = pmax(ifelse(is.na(WP_Group_UNC), 0, WP_Group_UNC),
ifelse(is.na(BSL_Group_UNC), 0, BSL_Group_UNC))
) |>
dplyr::arrange(Year_int, PAI_ID) |>
dplyr::select(-Year_int)
# -- Eqs.46/47/44/45/40/37/43: annual-rate calculation per PAI -------------
# Group by PAI_ID so the t0 stock lookup is PAI-scoped (each PAI's own baseline).
years_avail <- sort(unique(suppressWarnings(as.integer(delta$Year))))
t0_yr <- min(years_avail)
delta <- delta |>
dplyr::group_by(PAI_ID) |>
dplyr::mutate(
# Monitoring gap from the project baseline year (t0)
Monitoring_gap_yrs = suppressWarnings(as.integer(Year)) - t0_yr,
# Eq.46: annual SOC stock change -- BSL (Control zone) [t CO2e/yr]
# dC_BSL = (SOC_bsl,t - SOC_bsl,t0) / gap * A_PAI
dSOC_BSL_ha_yr = ifelse(Monitoring_gap_yrs > 0,
(BSL_SOC_ha - BSL_SOC_ha[Year == as.character(t0_yr)]) /
Monitoring_gap_yrs, 0),
dC_BSL_tCO2e_yr = dSOC_BSL_ha_yr * WP_area_ha,
# Eq.47: annual SOC stock change -- WP (Project zone) [t CO2e/yr]
# dC_WP = (SOC_wp,t - SOC_wp,t0) / gap * A_PAI
dSOC_WP_ha_yr = ifelse(Monitoring_gap_yrs > 0,
(WP_SOC_ha - WP_SOC_ha[Year == as.character(t0_yr)]) /
Monitoring_gap_yrs, 0),
dC_WP_tCO2e_yr = dSOC_WP_ha_yr * WP_area_ha,
# I_sign (Eqs.44/45): +1 if WP annual rate >= BSL annual rate, -1 otherwise
# When WP > BSL the UNC deduction reduces both values (conservative);
# when WP < BSL the sign reverses to add uncertainty to the loss.
I_sign = ifelse(dC_WP_tCO2e_yr >= dC_BSL_tCO2e_yr, 1L, -1L),
# Eq.44: UNC-adjusted annual BSL stock change
# dC_BSL_unc = dC_BSL * (1 - Group_UNC * I_sign)
dC_BSL_unc = dC_BSL_tCO2e_yr * (1 - Group_UNC_combined * I_sign),
# Eq.45: UNC-adjusted annual WP stock change
# dC_WP_unc = dC_WP * (1 - Group_UNC * I_sign)
dC_WP_unc = dC_WP_tCO2e_yr * (1 - Group_UNC_combined * I_sign),
# I_cumwp: switches ON (=1) once cumulative UNC-adjusted WP stock > 0
# (i.e., the project has accumulated a net carbon gain above baseline)
I_cumwp = as.integer(cumsum(dC_WP_unc) > 0),
# Eq.40: Gross carbon removals (t CO2e/yr)
# CR_t = I_cumwp x (MAX(0, dC_WP_unc) - MAX(0, dC_BSL_unc))
CR_t_tCO2e = I_cumwp *
(pmax(0, dC_WP_unc) - pmax(0, dC_BSL_unc)),
# Eq.37: Gross GHG emission reductions -- SOC component only [t CO2e/yr]
# Non-SOC components (fuel, lime, N2O) = 0; include via VCU model.
# ER_t = I_cumwp * (MIN(0,dC_WP_unc) - MIN(0,dC_BSL_unc))
# + (1-I_cumwp)* (MIN(0,dC_WP_unc) - MIN(0,dC_BSL_unc)
# + MAX(0,dC_WP_unc) - MAX(0,dC_BSL_unc))
ER_t_tCO2e = I_cumwp *
(pmin(0, dC_WP_unc) - pmin(0, dC_BSL_unc)) +
(1L - I_cumwp) *
(pmin(0, dC_WP_unc) - pmin(0, dC_BSL_unc) +
pmax(0, dC_WP_unc) - pmax(0, dC_BSL_unc)),
# Eq.43: Net ERR = ER_NET + CR_NET (leakage = 0; add via VCU model)
ERR_NET_tCO2e = ER_t_tCO2e + CR_t_tCO2e
) |>
dplyr::ungroup() |>
dplyr::arrange(suppressWarnings(as.integer(Year)), PAI_ID)
delta
} else data.frame()
cat("DELTA_NET_SOC rows:", nrow(delta_net_soc), "\n")## DELTA_NET_SOC rows: 12# ==============================================================================
# Step 7: Print group UNC summary
# ==============================================================================
cat("\n--- Group UNC Summary (applied to all PAIs) ---\n")##
## --- Group UNC Summary (applied to all PAIs) ---## Year Zone Group_Total_area_ha Group_AreaWt_SOC_tCO2e_ha Group_Pooled_Var
## 1 2022 Control 3218.101 218.7796 1.781098
## 2 2022 Project 3769.558 151.5702 20.246901
## 3 2024 Control 3218.101 225.0019 3.702000
## 4 2024 Project 3769.558 197.0401 12.057920
## n_PAIs Group_UNC_fraction Group_UNC_pct
## 1 6 0.002627314 0.2627314
## 2 6 0.012786167 1.2786167
## 3 6 0.003683043 0.3683043
## 4 6 0.007590257 0.7590257# ==============================================================================
# Styles
# ==============================================================================
wb <- openxlsx::createWorkbook()
hdr_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10, fontColour="white",
fgFill="#1F3864", textDecoration="bold", halign="center", wrapText=TRUE)
alt_sty <- openxlsx::createStyle(fgFill="#F2F6FC")
num4_sty <- openxlsx::createStyle(numFmt="0.0000", halign="right", fontName="Arial", fontSize=10)
num2_sty <- openxlsx::createStyle(numFmt="0.00", halign="right", fontName="Arial", fontSize=10)
num0_sty <- openxlsx::createStyle(numFmt="#,##0", halign="right", fontName="Arial", fontSize=10)
pct2_sty <- openxlsx::createStyle(numFmt="0.00%", halign="right", fontName="Arial", fontSize=10)
txt_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10, halign="left")
ttl_sty <- openxlsx::createStyle(fontName="Arial", fontSize=13,
fontColour="#1F3864", textDecoration="bold")
proj_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10,
fontColour="white", fgFill="#1F3864", textDecoration="bold")
grn_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10,
fontColour="white", fgFill="#375623", textDecoration="bold")
red_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10,
fontColour="white", fgFill="#C00000", textDecoration="bold")
grp_sty <- openxlsx::createStyle(fontName="Arial", fontSize=10, # group UNC highlight
fontColour="white", fgFill="#7030A0", textDecoration="bold")
write_ws <- function(wb, sname, df, title=NULL,
id_cols=3, freeze_col=4, num_fmt=NULL) {
dr <- if (!is.null(title)) {
openxlsx::writeData(wb, sname, title, startRow=1, startCol=1)
openxlsx::addStyle(wb, sname, ttl_sty, rows=1, cols=1)
3L
} else 2L
openxlsx::writeData(wb, sname, df, startRow=dr-1L, startCol=1,
headerStyle=hdr_sty, rowNames=FALSE,
borders="all", borderStyle="thin", borderColour="#CCCCCC")
nr <- nrow(df); nc <- ncol(df)
data_rows <- dr:(dr + nr - 1L)
even_rows <- data_rows[seq(2, nr, 2)]
if (length(even_rows))
openxlsx::addStyle(wb, sname, alt_sty, rows=even_rows, cols=1:nc,
gridExpand=TRUE, stack=TRUE)
if (id_cols > 0)
openxlsx::addStyle(wb, sname, txt_sty, rows=data_rows,
cols=1:min(id_cols, nc), gridExpand=TRUE, stack=TRUE)
num_idx <- which(sapply(df, is.numeric))
for (ci in num_idx) {
sty <- if (!is.null(num_fmt) && !is.null(num_fmt[[ci]])) num_fmt[[ci]] else num4_sty
openxlsx::addStyle(wb, sname, sty, rows=data_rows, cols=ci,
gridExpand=TRUE, stack=TRUE)
}
openxlsx::setColWidths(wb, sname, cols=1:nc,
widths=c(rep(16, min(id_cols,nc)), rep(14, max(0, nc-id_cols))))
openxlsx::freezePane(wb, sname, firstActiveRow=dr, firstActiveCol=freeze_col)
}
# Sheet order
snames <- c("README","EQ_REF","RUN_LOG","ALL_FD","ALL_ESM",
"SUMMARY_SOC","AREA_WT_SOC","DELTA_NET_SOC")
for (sn in snames) openxlsx::addWorksheet(wb, sn, gridLines=FALSE)
# ==============================================================================
# EQ_REF
# ==============================================================================
eq_ref <- data.frame(
Equation = c("Eq.46","Eq.47","Eq.74","Eq.44","Eq.45",
"Eq.37","Eq.38","Eq.39","Eq.40","Eq.41","Eq.42","Eq.43",
"Eq.75","Eq.76","Eq.77","Eq.78","Eq.79"),
Symbol = c("dC_BSL,s,t","dC_WP,s,t","UNC_t",
"dC_BSL_unc,t","dC_WP_unc,t",
"ER_t","ER_NET,t","LK_ER,t","CR_t","CR_NET,t","LK_CR,t","ERR_NET,t",
"BuER,t","BuCR,t","VCU_ER,t","VCU_CR,t","VCU_t"),
VM0042_Section = c("8.2","8.2","8.6","8.2","8.2",
"8.5","8.5","8.5","8.5","8.5","8.5","8.5",
"8.7","8.7","8.7","8.7","8.7"),
Formula = c(
"dC_BSL = (SOC_BSL_t - SOC_BSL_t-1) / verif_x",
"dC_WP = (SOC_WP_t - SOC_WP_t-1) / verif_x",
"UNC = min(1, sqrt(Var_aw) / |mean_SOC| * t_0.667) [applied at GROUP level in v1.8]",
"dC_BSL_unc = dC_BSL * (1 - Group_UNC * I_sign)",
"dC_WP_unc = dC_WP * (1 - Group_UNC * I_sign); I_sign=+1 if dC_WP>=dC_BSL else -1",
"ER = I_cumwp*(GHG+min(0,dC_WP_unc)-min(0,dC_BSL_unc)) + (1-I_cumwp)*same+max terms",
"ER_NET = ER - LK_ER",
"LK_ER = Leakage * ER / (ER+CR)",
"CR = I_cumwp * (max(0,dC_WP_unc) - max(0,dC_BSL_unc))",
"CR_NET = CR - LK_CR",
"LK_CR = Leakage * CR / (ER+CR)",
"ERR_NET = ER_NET + CR_NET",
"BuER = NPR_pct * [I*(min(0,WP)-min(0,BSL)) + (1-I)*all_terms]",
"BuCR = NPR_pct * I * (max(0,WP) - max(0,BSL))",
"VCU_ER = max(0, ER_NET - BuER)",
"VCU_CR = max(0, CR_NET - BuCR)",
"VCU = VCU_ER + VCU_CR"
),
Description = c(
"Annual SOC stock change -- BSL (Control) zone, t CO2e/ha/yr",
"Annual SOC stock change -- WP (Project) zone, t CO2e/ha/yr",
paste0("Uncertainty deduction fraction; pooled at GROUP (all-PAI) level in v1.8; t_0.667=", t0667),
"BSL SOC change adjusted for group uncertainty (conservative deduction)",
"WP SOC change adjusted for group uncertainty (conservative deduction)",
"Gross GHG emission reductions (SOC + non-SOC components)",
"Net GHG emission reductions after leakage deduction",
"Leakage allocated to reductions",
"Gross carbon removals (WP SOC increase above BSL)",
"Net carbon removals after leakage deduction",
"Leakage allocated to removals",
"Net ERR = net emission reductions + net carbon removals",
"Buffer pool deduction from reductions (non-permanence risk)",
"Buffer pool deduction from removals (non-permanence risk)",
"VCUs from reductions (after buffer)",
"VCUs from removals (after buffer)",
"Total annual VCU issuance"
),
R_Object = c("dSOC_BSL_ha_yr / dC_BSL_tCO2e_yr","dSOC_WP_ha_yr / dC_WP_tCO2e_yr",
"Group_UNC_fraction / Group_UNC_combined",
"dC_BSL_unc","dC_WP_unc",
"ER_t_tCO2e (SOC only)","-- see VCU model --",
"-- see VCU model --","CR_t_tCO2e","-- see VCU model --",
"-- see VCU model --","ERR_NET_tCO2e",
"-- see VCU model --","-- see VCU model --",
"-- see VCU model --","-- see VCU model --","-- see VCU model --"),
Where_Calculated = c(
"DELTA_NET_SOC","DELTA_NET_SOC","AREA_WT_SOC (Group_UNC_fraction col)",
"DELTA_NET_SOC","DELTA_NET_SOC",
"DELTA_NET_SOC (SOC only)","VCU Model","VCU Model",
"DELTA_NET_SOC","VCU Model","VCU Model","DELTA_NET_SOC / VCU Model",
"VCU Model","VCU Model","VCU Model","VCU Model","VCU Model"
),
stringsAsFactors = FALSE
)
openxlsx::writeData(wb, "EQ_REF",
paste0("VM0042 v2.2 Equation Reference | Project 4842 SAS | ",
"For IQC / VVB Audit | t0.667=", t0667, " | ACRES_TO_HA=", ACRES_TO_HA,
" | v1.8: UNC applied at GROUP level"),
startRow=1, startCol=1)
openxlsx::addStyle(wb, "EQ_REF", ttl_sty, rows=1, cols=1)
openxlsx::writeData(wb, "EQ_REF", eq_ref, startRow=3, startCol=1,
headerStyle=hdr_sty, rowNames=FALSE,
borders="all", borderStyle="thin", borderColour="#CCCCCC")
nr_eq <- nrow(eq_ref)
for (row_i in seq_len(nr_eq)) {
if (row_i %% 2 == 0)
openxlsx::addStyle(wb, "EQ_REF", alt_sty, rows=row_i+3, cols=1:ncol(eq_ref),
gridExpand=TRUE, stack=FALSE)
}
key_eq_rows <- which(eq_ref$Equation %in% c("Eq.46","Eq.47","Eq.74","Eq.43")) + 3L
openxlsx::addStyle(wb, "EQ_REF", proj_sty, rows=key_eq_rows, cols=1:ncol(eq_ref),
gridExpand=TRUE, stack=FALSE)
openxlsx::setColWidths(wb, "EQ_REF", cols=1:ncol(eq_ref), widths=c(10,18,14,80,55,25,25))
openxlsx::freezePane(wb, "EQ_REF", firstActiveRow=4, firstActiveCol=1)
# ==============================================================================
# README
# ==============================================================================
run_log_df <- if (length(run_log)>0) dplyr::bind_rows(run_log) else
data.frame(Status="none", stringsAsFactors=FALSE)
depth_sum <- if (length(sheet_depth_log)>0) dplyr::bind_rows(sheet_depth_log) else
data.frame()
depth_modes <- if (nrow(depth_sum)>0)
paste(unique(depth_sum$ESM_depths), collapse=" / ") else "n/a"
readme_df <- data.frame(
Field = c(
"Project","Script","Version","Method","Reference",
"Input file","Output file","Run date-time",
"PAIs","Years","Zones",
"","Sheets","README","EQ_REF","RUN_LOG","ALL_FD","ALL_ESM",
"SUMMARY_SOC","AREA_WT_SOC","DELTA_NET_SOC",
"","Assumptions",
"A1","A2","A3","A4","A5","A6","A7",
"","Column mapping (V1_0 file)",
"Sample ID","Rep_ID","Rep_old/Rep.1","Lower_cm(=0)","Upper_cm(=30)",
"","ESM depth detection",
"Mode","Detected depths",
"","Key VM0042 constants",
"t0.667","ACRES_TO_HA","G_CM2_TO_T_HA","C_TO_CO2e",
"","Changes vs v1.7",
"Change 1","Change 2","Change 3","Change 4","Change 5"
),
Value = c(
"4842 SAS -- Soil & Agriculture Sustainability",
"4842_ESM_SOC_ALL_PAIs_v1_8.Rmd","v1.8",
"Equivalent Soil Mass (ESM) -- Haden et al. 2020",
"von Haden AC, Yang WH, DeLucia EH (2020). Soil Biology & Biochemistry.",
input_file_name, basename(output_filename),
format(Sys.time(), "%Y-%m-%d %H:%M:%S %Z"),
paste(pai_sheets, collapse=", "),
paste(unique(run_log_df$Year[run_log_df$Status=="OK"]), collapse=", "),
"Project (WP) and Control (BSL) -- both processed",
"","Description",
"Version, method, assumptions, column mapping",
"VM0042 v2.2 equation reference for IQC/VVB audit",
"Processing status per PAI x Year x Zone (includes detected ESM_depths)",
"Combined fixed-depth data -- both zones",
"Combined ESM-corrected data -- both zones",
"Mean / n / SD (t CO2e/ha) + Mean_SOC_pct + SD_SOC_pct per PAI x Year x Zone x Strata",
paste0("Area-weighted SOC per PAI x Year x Zone; PAI_UNC and Group_UNC both shown;",
" Group_UNC used for SOC_net deduction"),
"Net SOC removal (CR) and reduction (ER) per PAI using Group_UNC_combined",
"","VM0042 v2.2 / VCS v4.7",
paste0("t0.667=",t0667," (one-sided t at 66.7% CI, large-n limit, VM0042 Eq.74)"),
"verif_period_yrs: auto-derived from (monitoring year - baseline year)",
"Strata_area in source file in ACRES; converted using NIST exact 1/2.47105381",
"Cum_SOC_g_cm2 * 100 = t C/ha (1 g/cm2=10000 g/m2/1e6 g/t*10000 m2/ha=100 t/ha)",
"t C/ha * (44/12) = t CO2e/ha (IPCC molar mass ratio)",
"Group_UNC_combined = max(WP_Group_UNC, BSL_Group_UNC) conservatively",
paste0("Group UNC pooled across all PAIs for same Year x Zone;",
" applied uniformly to every PAI (not individual PAI UNC)"),
"","Applied in prepare_columns()",
"-> ID (sample identifier)",
"-> Ref_ID (ESM reference profile for spline)",
"-> Rep (only renamed if Rep_old/Rep.1 present; V1_0 'Rep' kept as-is)",
"-> Upper_cm = 0 cm surface (file INVERTED; swap applied)",
"-> Lower_cm = 30 cm bottom (file INVERTED; swap applied)",
"","",
if (is.null(esm_depths_override)) "AUTO-DETECT (esm_depths_override = NULL)" else
paste("OVERRIDE:", paste(esm_depths_override, collapse=",")),
depth_modes,
"","",
as.character(t0667),
as.character(round(ACRES_TO_HA, 10)),
as.character(G_CM2_TO_T_HA),
as.character(round(C_TO_CO2e, 6)),
"","Relative to v1.5",
paste0("Auto-detect ESM depths per sheet from unique Lower_cm values;",
" works for 0-30 single layer and multi-increment profiles (0-15-30, etc.)"),
paste0("Group (project-wide) UNC pooled across all PAIs and applied uniformly;",
" per-PAI UNC retained for transparency in AREA_WT_SOC;",
" DELTA_NET_SOC uses Group_UNC_combined"),
paste0("FIX v1.8: Wrapped fd_out/esm_out in as.data.frame() to fix openxlsx tibble",
" write bug; ALL_FD and ALL_ESM now correctly populated"),
paste0("ADD v1.8: Mean_SOC_pct and SD_SOC_pct columns in SUMMARY_SOC tab",
" (ESM-corrected SOC percent per stratum per PAI per year)"),
paste0("FIX v1.8: DELTA_NET_SOC sorted by Year then Zone then PAI_ID",
" so all BSL/WP rows for a given year appear together")
),
stringsAsFactors = FALSE
)
openxlsx::writeData(wb, "README", readme_df, headerStyle=hdr_sty, rowNames=FALSE)
openxlsx::setColWidths(wb, "README", cols=1:2, widths=c(28, 80))
openxlsx::freezePane(wb, "README", firstActiveRow=2, firstActiveCol=1)
# ==============================================================================
# RUN_LOG
# ==============================================================================
openxlsx::writeData(wb, "RUN_LOG", run_log_df, headerStyle=hdr_sty, rowNames=FALSE)
openxlsx::setColWidths(wb, "RUN_LOG", cols=1:ncol(run_log_df),
widths=rep(15, ncol(run_log_df)))
openxlsx::freezePane(wb, "RUN_LOG", firstActiveRow=2, firstActiveCol=1)
# ==============================================================================
# ALL_FD and ALL_ESM
# ==============================================================================
all_fd <- if (length(all_fd_rows)>0) as.data.frame(dplyr::bind_rows(all_fd_rows)) else data.frame()
all_esm <- if (length(all_esm_rows)>0) as.data.frame(dplyr::bind_rows(all_esm_rows)) else data.frame()
if (nrow(all_fd)>0)
write_ws(wb, "ALL_FD", all_fd,
title="All PAIs & Years & Zones -- Fixed-Depth Data | 4842 SAS ESM v1.8",
id_cols=5, freeze_col=6)
if (nrow(all_esm)>0)
write_ws(wb, "ALL_ESM", all_esm,
title="All PAIs & Years & Zones -- ESM-Corrected Data | 4842 SAS ESM v1.8",
id_cols=5, freeze_col=6)
cat("ALL_FD rows:", nrow(all_fd), "| ALL_ESM rows:", nrow(all_esm), "\n")## ALL_FD rows: 933 | ALL_ESM rows: 906# ==============================================================================
# SUMMARY_SOC
# ==============================================================================
if (nrow(summary_soc) > 0) {
nc_s <- ncol(summary_soc)
nfmt_s <- vector("list", nc_s)
for (ci in which(sapply(summary_soc, is.numeric))) nfmt_s[[ci]] <- num2_sty
n_idx <- which(names(summary_soc) == "n_samples")
if (length(n_idx)) nfmt_s[[n_idx]] <- num0_sty
d_idx <- which(names(summary_soc) == "ESM_max_depth_cm")
if (length(d_idx)) nfmt_s[[d_idx]] <- num0_sty
pct_idx <- which(names(summary_soc) %in% c("Mean_SOC_pct","SD_SOC_pct"))
for (ci in pct_idx) nfmt_s[[ci]] <- num4_sty
write_ws(wb, "SUMMARY_SOC", summary_soc,
title=paste0("ESM SOC (t CO2e/ha) + SOC% | PAI x Year x Zone x Strata | ",
"VM0042 v2.2 Eq.46/47 inputs | Reps averaged per sample before statistics"),
id_cols=4, freeze_col=5, num_fmt=nfmt_s)
ctrl_rows <- which(summary_soc$Zone == "Control") + 2L
ctrl_sty <- openxlsx::createStyle(fgFill="#FFF2CC", fontName="Arial", fontSize=10)
if (length(ctrl_rows))
openxlsx::addStyle(wb, "SUMMARY_SOC", ctrl_sty, rows=ctrl_rows,
cols=1:nc_s, gridExpand=TRUE, stack=TRUE)
}
cat("SUMMARY_SOC rows:", nrow(summary_soc), "\n")## SUMMARY_SOC rows: 62# ==============================================================================
# AREA_WT_SOC
# Column highlight: Group_UNC columns in purple to distinguish from PAI_UNC
# ==============================================================================
if (nrow(area_wt_soc) > 0) {
nc_a <- ncol(area_wt_soc)
nfmt_a <- vector("list", nc_a)
for (ci in which(sapply(area_wt_soc, is.numeric))) nfmt_a[[ci]] <- num2_sty
unc_cols <- which(names(area_wt_soc) %in%
c("PAI_UNC_fraction","Group_UNC_fraction",
"Pooled_Var_PAI","Group_Pooled_Var"))
for (ci in unc_cols) nfmt_a[[ci]] <- num4_sty
pct_cols <- which(names(area_wt_soc) %in% c("PAI_UNC_pct","Group_UNC_pct"))
for (ci in pct_cols) nfmt_a[[ci]] <- num2_sty
n_ci <- which(names(area_wt_soc) == "n_samples_total")
if (length(n_ci)) nfmt_a[[n_ci]] <- num0_sty
write_ws(wb, "AREA_WT_SOC", area_wt_soc,
title=paste0("Area-Weighted SOC (t CO2e/ha) & Uncertainty | Both Zones | ",
"PAI_UNC = per-parcel (audit reference) | ",
"Group_UNC = project-wide pooled (used for deduction) | t0.667=", t0667),
id_cols=3, freeze_col=4, num_fmt=nfmt_a)
# Highlight GROUP UNC columns in purple
grp_unc_col_idx <- which(names(area_wt_soc) %in%
c("Group_UNC_fraction","Group_UNC_pct",
"Group_Total_area_ha","Group_AreaWt_SOC_tCO2e_ha",
"SOC_net_tCO2e_ha","SOC_net_tCO2e"))
if (length(grp_unc_col_idx)) {
all_data_rows <- (3):(2 + nrow(area_wt_soc))
openxlsx::addStyle(wb, "AREA_WT_SOC",
openxlsx::createStyle(fgFill="#EAD1FF", fontName="Arial", fontSize=10),
rows=all_data_rows, cols=grp_unc_col_idx,
gridExpand=TRUE, stack=TRUE)
# Header row highlight for group cols (row 2 = header when title in row 1)
openxlsx::addStyle(wb, "AREA_WT_SOC",
openxlsx::createStyle(fgFill="#7030A0", fontColour="white",
fontName="Arial", fontSize=10, textDecoration="bold",
halign="center", wrapText=TRUE),
rows=2, cols=grp_unc_col_idx, gridExpand=TRUE, stack=FALSE)
}
# PROJECT rows dark blue
proj_rows_a <- which(area_wt_soc$PAI_ID == "PROJECT") + 2L
if (length(proj_rows_a))
openxlsx::addStyle(wb, "AREA_WT_SOC", proj_sty, rows=proj_rows_a,
cols=1:nc_a, gridExpand=TRUE, stack=FALSE)
# Control rows yellow
ctrl_rows_a <- which(area_wt_soc$Zone == "Control" &
area_wt_soc$PAI_ID != "PROJECT") + 2L
ctrl_sty2 <- openxlsx::createStyle(fgFill="#FFF2CC", fontName="Arial", fontSize=10)
if (length(ctrl_rows_a))
openxlsx::addStyle(wb, "AREA_WT_SOC", ctrl_sty2, rows=ctrl_rows_a,
cols=1:nc_a, gridExpand=TRUE, stack=TRUE)
}
cat("AREA_WT_SOC rows:", nrow(area_wt_soc), "\n")## AREA_WT_SOC rows: 28# ==============================================================================
# DELTA_NET_SOC
# ==============================================================================
if (nrow(delta_net_soc) > 0) {
nc_d <- ncol(delta_net_soc)
nfmt_d <- vector("list", nc_d)
for (ci in which(sapply(delta_net_soc, is.numeric))) nfmt_d[[ci]] <- num2_sty
# 4-dp for UNC fractions
unc_d <- which(names(delta_net_soc) %in%
c("WP_Group_UNC","BSL_Group_UNC","Group_UNC_combined"))
if (length(unc_d)) for (ci in unc_d) nfmt_d[[ci]] <- num4_sty
# Integer for gap and indicators
int_d <- which(names(delta_net_soc) %in% c("Monitoring_gap_yrs","I_sign","I_cumwp"))
if (length(int_d)) for (ci in int_d) nfmt_d[[ci]] <- num0_sty
write_ws(wb, "DELTA_NET_SOC", delta_net_soc,
title=paste0(
"Net SOC Removal (CR_t) and Reduction (ER_t) per PAI | ",
"VM0042 v2.2 Eqs.37,40,43,44,45,46,47,74 | ",
"Annual rates (Eq.46/47) -> UNC-adjusted (Eq.44/45) -> CR/ER/ERR_NET (Eq.40/37/43) | ",
"Non-SOC GHG (fuel, lime, N2O) = 0 -- add via VCU Model"),
id_cols=2, freeze_col=3, num_fmt=nfmt_d)
# Green rows where CR>0 (net removal), red where ER>0 (reduction)
if ("CR_t_tCO2e" %in% names(delta_net_soc)) {
cr_rows <- which(!is.na(delta_net_soc$CR_t_tCO2e) &
delta_net_soc$CR_t_tCO2e > 0) + 2L
er_rows <- which(!is.na(delta_net_soc$ER_t_tCO2e) &
delta_net_soc$ER_t_tCO2e > 0) + 2L
if (length(cr_rows))
openxlsx::addStyle(wb, "DELTA_NET_SOC", grn_sty, rows=cr_rows,
cols=1:nc_d, gridExpand=TRUE, stack=FALSE)
if (length(er_rows))
openxlsx::addStyle(wb, "DELTA_NET_SOC", red_sty, rows=er_rows,
cols=1:nc_d, gridExpand=TRUE, stack=FALSE)
}
}
cat("DELTA_NET_SOC rows:", nrow(delta_net_soc), "\n")## DELTA_NET_SOC rows: 12# ==============================================================================
# SAVE
# ==============================================================================
tryCatch({
openxlsx::saveWorkbook(wb, output_filename, overwrite=TRUE)
cat("\n========================================================\n")
cat("SAVED:", output_filename, "\n")
cat("Sheets:", paste(openxlsx::sheets(wb), collapse=" | "), "\n")
cat("========================================================\n")
}, error = function(e) warning("Save failed: ", conditionMessage(e)))##
## ========================================================
## SAVED: C:/Users/Quan/Carbon Friendly/Carbon Friendly - Verra/Project 4842 Sentry Ag Services Carbon Project/Validation/DATA/ESM_SOC_ALL_PAIs_20260402_185811.xlsx
## Sheets: README | EQ_REF | RUN_LOG | ALL_FD | ALL_ESM | SUMMARY_SOC | AREA_WT_SOC | DELTA_NET_SOC
## ========================================================## ### Group (Project-Wide) Uncertainty Applied to All PAIs| Year | Zone | Group_Total_area_ha | Group_AreaWt_SOC_tCO2e_ha | Group_Pooled_Var | n_PAIs | Group_UNC_fraction | Group_UNC_pct |
|---|---|---|---|---|---|---|---|
| 2022 | Control | 3218.101 | 218.7796 | 1.7811 | 6 | 0.0026 | 0.2627 |
| 2022 | Project | 3769.558 | 151.5702 | 20.2469 | 6 | 0.0128 | 1.2786 |
| 2024 | Control | 3218.101 | 225.0019 | 3.7020 | 6 | 0.0037 | 0.3683 |
| 2024 | Project | 3769.558 | 197.0401 | 12.0579 | 6 | 0.0076 | 0.7590 |
| PAI_ID | Year | Zone | Strata | Strata_area_ha | ESM_max_depth_cm | n_samples | Mean_SOC_tCO2e_ha | SD_SOC_tCO2e_ha | SE_SOC_tCO2e_ha | Mean_SOC_tC_ha | SD_SOC_tC_ha | Mean_SOC_pct | SD_SOC_pct |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| PAI11 | 2022 | Control | Stratum 0 | 132.0 | 30 | 1 | 291.12 | NA | NA | 79.39511 | NA | 2.0100000 | NA |
| PAI11 | 2022 | Control | Stratum 2 | 148.9 | 30 | 1 | 291.12 | NA | NA | 79.39502 | NA | 2.0100000 | NA |
| PAI11 | 2022 | Project | Stratum 0 | 132.0 | 30 | 5 | 189.86 | 70.90 | 31.71 | 51.78031 | 19.337486 | 1.2520000 | 0.5333104 |
| PAI11 | 2022 | Project | Stratum 1 | 100.2 | 30 | 3 | 258.45 | 50.83 | 29.35 | 70.48600 | 13.862739 | 1.7700000 | 0.3296968 |
| PAI11 | 2022 | Project | Stratum 2 | 148.9 | 30 | 4 | 188.98 | 7.56 | 3.78 | 51.53938 | 2.062870 | 1.2650000 | 0.0900000 |
| PAI11 | 2024 | Control | Stratum 0 | 132.0 | 30 | 1 | 282.43 | NA | NA | 77.02511 | NA | 1.9500000 | NA |
| PAI11 | 2024 | Control | Stratum 2 | 148.9 | 30 | 1 | 282.43 | NA | NA | 77.02502 | NA | 1.9500000 | NA |
| PAI11 | 2024 | Project | Stratum 0 | 132.0 | 30 | 5 | 219.09 | 90.44 | 40.44 | 59.75051 | 24.664455 | 1.4560000 | 0.6577462 |
| PAI11 | 2024 | Project | Stratum 1 | 100.2 | 30 | 3 | 288.06 | 50.31 | 29.05 | 78.56134 | 13.720416 | 1.9733333 | 0.3257811 |
| PAI11 | 2024 | Project | Stratum 2 | 148.9 | 30 | 4 | 230.03 | 13.62 | 6.81 | 62.73463 | 3.713413 | 1.5400000 | 0.1283225 |
| PAI15 | 2022 | Control | Stratum 1 | 858.4 | 30 | 1 | 247.75 | NA | NA | 67.56798 | NA | 1.6400000 | NA |
| PAI15 | 2022 | Control | Stratum 2 | 1101.7 | 30 | 1 | 247.75 | NA | NA | 67.56806 | NA | 1.6400000 | NA |
| PAI15 | 2022 | Project | Stratum 0 | 83.3 | 30 | 3 | 168.85 | 68.84 | 39.74 | 46.04900 | 18.774089 | 1.0500000 | 0.4396590 |
| PAI15 | 2022 | Project | Stratum 1 | 858.4 | 30 | 20 | 152.02 | 48.46 | 10.84 | 41.45955 | 13.216482 | 0.9535000 | 0.3328391 |
| PAI15 | 2022 | Project | Stratum 2 | 1101.7 | 30 | 27 | 125.71 | 50.11 | 9.64 | 34.28422 | 13.665674 | 0.7792593 | 0.3578668 |
| PAI15 | 2024 | Control | Stratum 1 | 858.4 | 30 | 1 | 238.69 | NA | NA | 65.09598 | NA | 1.5800000 | NA |
| PAI15 | 2024 | Control | Stratum 2 | 1101.7 | 30 | 1 | 238.69 | NA | NA | 65.09606 | NA | 1.5800000 | NA |
| PAI15 | 2024 | Project | Stratum 0 | 83.3 | 30 | 7 | 228.37 | 45.90 | 17.35 | 62.28329 | 12.517007 | 1.4671429 | 0.3015358 |
| PAI15 | 2024 | Project | Stratum 1 | 858.4 | 30 | 32 | 200.14 | 43.80 | 7.74 | 54.58478 | 11.946430 | 1.2721875 | 0.3090228 |
| PAI15 | 2024 | Project | Stratum 2 | 1101.7 | 30 | 34 | 176.65 | 40.08 | 6.87 | 48.17753 | 10.930268 | 1.1023529 | 0.2939600 |
| PAI16 | 2022 | Control | Stratum 0 | NA | 30 | 1 | 253.02 | NA | NA | 69.00546 | NA | 1.5700000 | NA |
| PAI16 | 2022 | Control | Stratum 1 | 38.4 | 30 | 1 | 253.02 | NA | NA | 69.00553 | NA | 1.5700000 | NA |
| PAI16 | 2022 | Control | Stratum 2 | 58.5 | 30 | 1 | 253.02 | NA | NA | 69.00556 | NA | 1.5700000 | NA |
| PAI16 | 2022 | Project | Stratum 1 | 38.4 | 30 | 2 | 223.62 | 15.37 | 10.87 | 60.98764 | 4.192687 | 1.4400000 | 0.0989949 |
| PAI16 | 2022 | Project | Stratum 2 | 58.5 | 30 | 2 | 223.62 | 15.37 | 10.87 | 60.98764 | 4.192687 | 1.4400000 | 0.0989949 |
| PAI16 | 2024 | Control | Stratum 0 | NA | 30 | 1 | 251.41 | NA | NA | 68.56593 | NA | 1.5600000 | NA |
| PAI16 | 2024 | Control | Stratum 1 | 38.4 | 30 | 1 | 251.41 | NA | NA | 68.56601 | NA | 1.5600000 | NA |
| PAI16 | 2024 | Control | Stratum 2 | 58.5 | 30 | 1 | 251.41 | NA | NA | 68.56603 | NA | 1.5600000 | NA |
| PAI16 | 2024 | Project | Stratum 1 | 38.4 | 30 | 2 | 243.03 | 34.04 | 24.07 | 66.28171 | 9.283806 | 1.5650000 | 0.2192031 |
| PAI16 | 2024 | Project | Stratum 2 | 58.5 | 30 | 2 | 243.03 | 34.04 | 24.07 | 66.28171 | 9.283806 | 1.5650000 | 0.2192031 |
| PAI3 | 2022 | Control | Stratum 0 | 209.2 | 30 | 1 | 126.23 | NA | NA | 34.42501 | NA | 0.8100000 | NA |
| PAI3 | 2022 | Control | Stratum 1 | 328.0 | 30 | 2 | 139.32 | 18.52 | 13.10 | 37.99650 | 5.050851 | 0.8850000 | 0.1060660 |
| PAI3 | 2022 | Project | Stratum 0 | 209.2 | 30 | 5 | 161.75 | 23.26 | 10.40 | 44.11440 | 6.342769 | 1.0320000 | 0.1583351 |
| PAI3 | 2022 | Project | Stratum 1 | 328.0 | 30 | 6 | 138.13 | 13.21 | 5.39 | 37.67050 | 3.602834 | 0.8783333 | 0.0870441 |
| PAI3 | 2022 | Project | Stratum 2 | 270.3 | 30 | 7 | 176.05 | 26.33 | 9.95 | 48.01457 | 7.181778 | 1.1485714 | 0.1546424 |
| PAI3 | 2024 | Control | Stratum 0 | 209.2 | 30 | 1 | 187.00 | NA | NA | 51.00001 | NA | 1.2000000 | NA |
| PAI3 | 2024 | Control | Stratum 1 | 328.0 | 30 | 2 | 168.12 | 26.70 | 18.88 | 45.85100 | 7.281801 | 1.0700000 | 0.1838478 |
| PAI3 | 2024 | Project | Stratum 0 | 209.2 | 30 | 5 | 186.54 | 30.47 | 13.63 | 50.87320 | 8.310263 | 1.1900000 | 0.2040833 |
| PAI3 | 2024 | Project | Stratum 1 | 328.0 | 30 | 6 | 161.97 | 8.82 | 3.60 | 44.17233 | 2.404995 | 1.0300000 | 0.0622896 |
| PAI3 | 2024 | Project | Stratum 2 | 270.3 | 30 | 7 | 201.49 | 23.44 | 8.86 | 54.95242 | 6.392877 | 1.3157143 | 0.1387873 |
| PAI5 | 2022 | Control | Stratum 0 | 87.6 | 30 | 1 | 127.16 | NA | NA | 34.67995 | NA | 0.8000000 | NA |
| PAI5 | 2022 | Control | Stratum 1 | 56.3 | 30 | 1 | 131.74 | NA | NA | 35.92777 | NA | 0.8600000 | NA |
| PAI5 | 2022 | Control | Stratum 2 | 105.3 | 30 | 1 | 125.14 | NA | NA | 34.12938 | NA | 0.8600000 | NA |
| PAI5 | 2022 | Project | Stratum 0 | 87.6 | 30 | 2 | 115.72 | 5.29 | 3.74 | 31.56091 | 1.441884 | 0.7250000 | 0.1060660 |
| PAI5 | 2022 | Project | Stratum 1 | 56.3 | 30 | 2 | 92.89 | 27.00 | 19.09 | 25.33384 | 7.364523 | 0.5600000 | 0.1272792 |
| PAI5 | 2022 | Project | Stratum 2 | 105.3 | 30 | 3 | 135.03 | 13.65 | 7.88 | 36.82511 | 3.723323 | 0.8533333 | 0.0472582 |
| PAI5 | 2024 | Control | Stratum 0 | 87.6 | 30 | 1 | 170.08 | NA | NA | 46.38444 | NA | 1.0700000 | NA |
| PAI5 | 2024 | Control | Stratum 1 | 56.3 | 30 | 1 | 194.54 | NA | NA | 53.05614 | NA | 1.2700000 | NA |
| PAI5 | 2024 | Control | Stratum 2 | 105.3 | 30 | 1 | 184.80 | NA | NA | 50.40036 | NA | 1.2700000 | NA |
| PAI5 | 2024 | Project | Stratum 0 | 87.6 | 30 | 1 | 238.92 | NA | NA | 65.16096 | NA | 1.6000000 | NA |
| PAI5 | 2024 | Project | Stratum 1 | 56.3 | 30 | 1 | 204.12 | NA | NA | 55.66860 | NA | 1.3000000 | NA |
| PAI5 | 2024 | Project | Stratum 2 | 105.3 | 30 | 3 | 215.65 | 24.57 | 14.18 | 58.81456 | 6.699599 | 1.3733333 | 0.2300725 |
| PAI9 | 2022 | Control | Stratum 0 | 79.8 | 30 | 1 | 88.80 | NA | NA | 24.21884 | NA | 0.6100000 | NA |
| PAI9 | 2022 | Control | Stratum 2 | 14.1 | 30 | 1 | 88.80 | NA | NA | 24.21885 | NA | 0.6100000 | NA |
| PAI9 | 2022 | Project | Stratum 0 | 79.8 | 30 | 2 | 171.99 | 88.61 | 62.66 | 46.90727 | 24.165985 | 1.1750000 | 0.7141778 |
| PAI9 | 2022 | Project | Stratum 1 | 97.6 | 30 | 4 | 162.06 | 58.89 | 29.45 | 44.19802 | 16.061732 | 1.1150000 | 0.4670831 |
| PAI9 | 2022 | Project | Stratum 2 | 14.1 | 30 | 1 | 119.27 | NA | NA | 32.52823 | NA | 0.8000000 | NA |
| PAI9 | 2024 | Control | Stratum 0 | 79.8 | 30 | 1 | 138.30 | NA | NA | 37.71787 | NA | 0.9500000 | NA |
| PAI9 | 2024 | Control | Stratum 2 | 14.1 | 30 | 1 | 138.30 | NA | NA | 37.71789 | NA | 0.9500000 | NA |
| PAI9 | 2024 | Project | Stratum 0 | 79.8 | 30 | 2 | 189.97 | 7.87 | 5.57 | 51.81082 | 2.147498 | 1.2650000 | 0.1909188 |
| PAI9 | 2024 | Project | Stratum 1 | 97.6 | 30 | 4 | 225.11 | 62.20 | 31.10 | 61.39231 | 16.963625 | 1.5350000 | 0.4912908 |
| PAI9 | 2024 | Project | Stratum 2 | 14.1 | 30 | 1 | 202.76 | NA | NA | 55.29799 | NA | 1.3600000 | NA |
| PAI_ID | Year | Zone | ESM_max_depth_cm | Total_area_ha | AreaWt_SOC_tCO2e_ha | Pooled_Var_PAI | n_samples_total | n_strata | n_strata_no_SD | PAI_UNC_fraction | PAI_UNC_pct | Total_SOC_tCO2e | PAI_Notes | Group_UNC_fraction | Group_UNC_pct | Group_Pooled_Var | Group_Total_area_ha | Group_AreaWt_SOC_tCO2e_ha | SOC_net_tCO2e_ha | SOC_net_tCO2e |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| PAI11 | 2022 | Control | 30 | 280.9044 | 291.1152 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 81775.56 | 2 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 290.3504 | 81560.712 |
| PAI15 | 2022 | Control | 30 | 1960.0504 | 247.7494 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 485601.38 | 2 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 247.0985 | 484325.552 |
| PAI16 | 2022 | Control | 30 | 96.8372 | 253.0203 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 24501.79 | 2 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 252.3556 | 24437.415 |
| PAI3 | 2022 | Control | 30 | 537.1838 | 134.2201 | 63.9207 | 3 | 2 | 1 | 0.0257 | 2.5655 | 72100.84 | 1 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 133.8674 | 71911.412 |
| PAI5 | 2022 | Control | 30 | 249.2418 | 127.3394 | 0.0000 | 3 | 3 | 3 | 0.0000 | 0.0000 | 31738.31 | 3 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 127.0048 | 31654.920 |
| PAI9 | 2022 | Control | 30 | 93.8830 | 88.8024 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 8337.04 | 2 strata n=1 (SD=NA) | 0.0026 | 0.2627 | 1.7811 | 3218.101 | 218.7796 | 88.5691 | 8315.136 |
| PROJECT | 2022 | Control | 30 | 3218.1007 | 218.7796 | 1.7811 | 14 | 13 | 12 | 0.0026 | 0.2627 | 704054.92 | PROJECT row: group pooled | 0.0026 | 0.2627 | NA | 3218.101 | 218.7796 | 218.2048 | 702205.147 |
| PAI11 | 2022 | Project | 30 | 381.1532 | 207.5555 | 182.3171 | 12 | 3 | 0 | 0.0280 | 2.8019 | 79110.42 | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 204.9016 | 78098.906 | |
| PAI15 | 2022 | Project | 30 | 2043.3347 | 138.5193 | 50.3770 | 50 | 3 | 0 | 0.0221 | 2.2069 | 283041.29 | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 136.7482 | 279422.277 | |
| PAI16 | 2022 | Project | 30 | 96.8372 | 223.6214 | 61.6212 | 4 | 2 | 0 | 0.0151 | 1.5119 | 21654.87 | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 220.7621 | 21377.989 | |
| PAI3 | 2022 | Project | 30 | 807.5057 | 156.9439 | 23.1612 | 18 | 3 | 0 | 0.0132 | 1.3207 | 126733.12 | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 154.9372 | 125112.693 | |
| PAI5 | 2022 | Project | 30 | 249.2418 | 118.7273 | 31.4047 | 7 | 3 | 0 | 0.0203 | 2.0329 | 29591.81 | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 117.2092 | 29213.444 | |
| PAI9 | 2022 | Project | 30 | 191.4851 | 163.0467 | 906.7311 | 7 | 3 | 1 | 0.0795 | 7.9543 | 31221.02 | 1 strata n=1 (SD=NA) | 0.0128 | 1.2786 | 20.2469 | 3769.558 | 151.5702 | 160.9620 | 30821.822 |
| PROJECT | 2022 | Project | 30 | 3769.5577 | 151.5702 | 20.2469 | 98 | 17 | 1 | 0.0128 | 1.2786 | 571352.54 | PROJECT row: group pooled | 0.0128 | 1.2786 | NA | 3769.558 | 151.5702 | 149.6322 | 564047.131 |
| PAI11 | 2024 | Control | 30 | 280.9044 | 282.4252 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 79334.50 | 2 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 281.3850 | 79042.307 |
| PAI15 | 2024 | Control | 30 | 1960.0504 | 238.6854 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 467835.46 | 2 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 237.8063 | 466112.405 |
| PAI16 | 2024 | Control | 30 | 96.8372 | 251.4087 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 24345.73 | 2 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 250.4828 | 24256.060 |
| PAI3 | 2024 | Control | 30 | 537.1838 | 175.4736 | 132.8587 | 3 | 2 | 1 | 0.0283 | 2.8292 | 94261.58 | 1 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 174.8273 | 93914.406 |
| PAI5 | 2024 | Control | 30 | 249.2418 | 181.8221 | 0.0000 | 3 | 3 | 3 | 0.0000 | 0.0000 | 45317.68 | 3 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 181.1525 | 45150.774 |
| PAI9 | 2024 | Control | 30 | 93.8830 | 138.2989 | 0.0000 | 2 | 2 | 2 | 0.0000 | 0.0000 | 12983.92 | 2 strata n=1 (SD=NA) | 0.0037 | 0.3683 | 3.7020 | 3218.101 | 225.0019 | 137.7895 | 12936.095 |
| PROJECT | 2024 | Control | 30 | 3218.1007 | 225.0019 | 3.7020 | 14 | 13 | 12 | 0.0037 | 0.3683 | 724078.86 | PROJECT row: group pooled | 0.0037 | 0.3683 | NA | 3218.101 | 225.0019 | 224.1732 | 721412.047 |
| PAI11 | 2024 | Project | 30 | 381.1532 | 241.5013 | 261.5610 | 12 | 3 | 0 | 0.0288 | 2.8843 | 92048.97 | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 239.6682 | 91350.298 | |
| PAI15 | 2024 | Project | 30 | 2043.3347 | 188.6282 | 24.8140 | 73 | 3 | 0 | 0.0114 | 1.1374 | 385430.56 | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 187.1965 | 382505.046 | |
| PAI16 | 2024 | Project | 30 | 96.8372 | 243.0329 | 302.1325 | 4 | 2 | 0 | 0.0308 | 3.0804 | 23534.64 | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 241.1883 | 23356.002 | |
| PAI3 | 2024 | Project | 30 | 807.5057 | 181.5633 | 23.4003 | 18 | 3 | 0 | 0.0115 | 1.1475 | 146613.41 | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 180.1852 | 145500.580 | |
| PAI5 | 2024 | Project | 30 | 249.2418 | 221.2313 | 35.9331 | 5 | 3 | 2 | 0.0117 | 1.1670 | 55140.11 | 2 strata n=1 (SD=NA) | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 219.5521 | 54721.581 |
| PAI9 | 2024 | Project | 30 | 191.4851 | 208.8220 | 256.6674 | 7 | 3 | 1 | 0.0330 | 3.3043 | 39986.30 | 1 strata n=1 (SD=NA) | 0.0076 | 0.7590 | 12.0579 | 3769.558 | 197.0401 | 207.2370 | 39682.795 |
| PROJECT | 2024 | Project | 30 | 3769.5577 | 197.0401 | 12.0579 | 119 | 17 | 3 | 0.0076 | 0.7590 | 742754.00 | PROJECT row: group pooled | 0.0076 | 0.7590 | NA | 3769.558 | 197.0401 | 195.5445 | 737116.302 |
| PAI_ID | Year | WP_area_ha | WP_SOC_ha | BSL_SOC_ha | Delta_SOC_ha | Net_SOC_tCO2e | Group_UNC_combined | Monitoring_gap_yrs | dC_WP_tCO2e_yr | dC_BSL_tCO2e_yr | I_sign | dC_WP_unc | dC_BSL_unc | I_cumwp | CR_t_tCO2e | ER_t_tCO2e | ERR_NET_tCO2e |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| PAI11 | 2022 | 381.15 | 207.56 | 291.12 | -83.56 | -31849.07 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI15 | 2022 | 2043.33 | 138.52 | 247.75 | -109.23 | -223193.73 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI16 | 2022 | 96.84 | 223.62 | 253.02 | -29.40 | -2846.92 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI3 | 2022 | 807.51 | 156.94 | 134.22 | 22.72 | 18349.65 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI5 | 2022 | 249.24 | 118.73 | 127.34 | -8.61 | -2146.50 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI9 | 2022 | 191.49 | 163.05 | 88.80 | 74.24 | 14216.68 | 0.01 | 0 | 0.00 | 0.00 | 1 | 0.00 | 0.00 | 0 | 0.00 | 0.00 | 0.00 |
| PAI11 | 2024 | 381.15 | 241.50 | 282.43 | -40.92 | -15598.29 | 0.01 | 2 | 6469.27 | -1656.11 | 1 | 6420.17 | -1643.54 | 1 | 6420.17 | 1643.54 | 8063.71 |
| PAI15 | 2024 | 2043.33 | 188.63 | 238.69 | -50.06 | -102283.65 | 0.01 | 2 | 51194.64 | -9260.40 | 1 | 50806.06 | -9190.11 | 1 | 50806.06 | 9190.11 | 59996.17 |
| PAI16 | 2024 | 96.84 | 243.03 | 251.41 | -8.38 | -811.09 | 0.01 | 2 | 939.88 | -78.03 | 1 | 932.75 | -77.44 | 1 | 932.75 | 77.44 | 1010.19 |
| PAI3 | 2024 | 807.51 | 181.56 | 175.47 | 6.09 | 4917.48 | 0.01 | 2 | 9940.14 | 16656.23 | -1 | 10015.59 | 16782.66 | 1 | -6767.07 | 0.00 | -6767.07 |
| PAI5 | 2024 | 249.24 | 221.23 | 181.82 | 39.41 | 9822.43 | 0.01 | 2 | 12774.15 | 6789.69 | 1 | 12677.19 | 6738.15 | 1 | 5939.04 | 0.00 | 5939.04 |
| PAI9 | 2024 | 191.49 | 208.82 | 138.30 | 70.52 | 13504.13 | 0.01 | 2 | 4382.64 | 4738.92 | -1 | 4415.91 | 4774.89 | 1 | -358.98 | 0.00 | -358.98 |
## Script : 4842_ESM_SOC_ALL_PAIs_v1_8.Rmd## Version : v1.8## R : R version 4.4.1 (2024-06-14 ucrt)## Run : 2026-04-02 18:58:24 AEST## Input : C:/Users/Quan/Carbon Friendly/Carbon Friendly - Verra/Project 4842 Sentry Ag Services Carbon Project/Validation/DATA/4842_SOC_data_for_ESM_V1_0.xlsx## Output : C:/Users/Quan/Carbon Friendly/Carbon Friendly - Verra/Project 4842 Sentry Ag Services Carbon Project/Validation/DATA/ESM_SOC_ALL_PAIs_20260402_185811.xlsx## PAIs : PAI3, PAI5, PAI9, PAI11, PAI15, PAI16## t0.667 : 0.4307 (VM0042 Eq.74)## ESM depths : AUTO-DETECT per sheet## UNC mode : GROUP (project-wide pooled, applied uniformly to all PAIs)## Unit : g/cm2 * 100 = t C/ha; * 3.6667 = t CO2e/ha