Chapter 6 Comparison with VRI R1
The comparison process was designed to assess the consistency between the Beta version of the 2026 Provincial NSYT dataset and the currently published VRI R1. The workflow begins by loading the provincial NSYT yield table and joining it to the R1 inventory database using a common feature_id. Species, age, and volume information from both sources are combined to create a common analytical framework. To facilitate a direct comparison, stand ages are standardized, volume estimates are derived from the R1 inventory, and VDYP yield projections are interpolated to match the specific ages represented in the inventory. This allows both datasets to be evaluated using comparable volume estimates for the same forest features.
A difference metric is calculated to quantify the variation between the inventory-derived volume (volume from VRI R1 layer) and the NSYT volume estimate (interpolated NSYT volume). The resulting dataset provides a feature-by-feature comparison of the two yield systems and supports statistical evaluation of differences across the provincial land base. This process serves as a validation and quality assurance exercise, helping identify areas where modelling assumptions, age treatments, interpolation methods, or source inventory data may contribute to differences between the existing R1 inventory volumes and the proposed provincial NSYT yield dataset.
The following contrasts the Beta version of the 2026 NSYT dataset with the currently published production year 2025 (projection year Dec 2024) veg_comp_lyr_r1_poly_2024 (R1). The R1 contains 6872386 features, 90% of which (6209645) has a species identified and (spec_cd_1) and 89% (6140967) remain after initial cleaning.
library(DBI)
library(RPostgreSQL)
library(tidyverse)
library(terra)
library(R.utils)
library(tidytext)
library(arsenal)
library(dlookr)
library(gdalUtilities)
library(glue)
library(janitor)
library(knitr)
library(kableExtra)
library(tidytext)
library(fs)
library(data.table)
library(glue)
library(janitor)
library(dplyr)
rm(list = ls())
options(scipen = 999)
set.seed(123)
# source("log_files.R")
knitr::opts_chunk$set(
echo = TRUE, # Show code in the output
fig.width = 4, # Set the width of the figures (in inches)
fig.height = 4, # Set the height of the figures (in inches)
collapse = TRUE, # Collapsing code and output
# out.width = '250px', # Standardize output image width
# out.height = '250px', # Standardize output image height
message = FALSE, # Hide messages
warning = FALSE # Hide warnings
)6.1 Get data
The following snippet loads 2026 provincial NSYT from an RDS file containing 6140967 features, and joins it with the data pulled from the VRI vegetation composition R1 layer based on feature_id. Various new columns are created, including volumes from VRI R1 layer(vri_vol), capped age (vri_age250), and interpolated volumes (vdyp_vol) based on age ranges.
schema <- 'msyt_2025'
opt <- paste0("-c search_path=",schema)
user_name <- 'results'
database <- 'msyt'
db<- dbConnect('PostgreSQL',dbname=database,user=user_name,options = opt)
prov<-readRDS("./yld_trim.rds")
names(prov) <- tolower(names(prov))
sql<-glue_sql("select
feature_id,
proj_age_1 as age,
species_cd_1 as leading_species,
species_cd_2 as second_species,
live_stand_volume_125 as r1_vol125,
live_stand_volume_175 as r1_vol175,
dead_stand_volume_125 as r1_dvol125
from veg_comp_fdw_full",.con = db)
r1<-dbGetQuery(db,sql)
dbDisconnect(db)
##############################################################################
# #
prov_df<-inner_join(prov,r1,join_by(feature_id))
names(prov_df) <- tolower(names(prov_df))
prov_df2 <- prov_df %>%
mutate(
vri_vol = case_when(
grepl("^[PAEDMW]", leading_species) ~ r1_vol125,
leading_species %in% c("QG", "ZH", "VB", "RA", "XH") ~ r1_vol125,
TRUE ~ r1_vol175 # Otherwise, assign r1_vol175
),
across(where(is.numeric), ~replace_na(., 0)),
vri_age250 = pmin(age, 250),
age_low = floor(vri_age250/10) * 10,
age_high = ceiling(vri_age250/10) * 10,
p_low = (age_high - vri_age250) / 10
) %>%
group_by(feature_id) %>%
mutate(
l_vol = prj_vol_dwb[match(age_low, prj_total_age)],
u_vol = prj_vol_dwb[match(age_high, prj_total_age)],
vdyp_vol = (p_low * l_vol) + ((1 - p_low) * u_vol)
) %>%
select(feature_id, age, vri_vol,l_vol, u_vol, vdyp_vol)%>%
slice(1)%>%
ungroup()6.2 Interpolated Volume Comparison
The projected volume is linearly interpolated to the exact VRI age using the surrounding decadal projection ages.
VRI Volume: vri_vol is VRI volume (live_stand_volume_125 or live_stand_volume_175) selected based on the leading species. Volume Interpolation: vdyp_vol is calculated by linearly interpolating between the projection volumes at the lower and upper decadal ages. Volume Difference: vol_dif is calculated as vri_vol - vdyp_vol.
| described_variables | n | na | mean | sd | se_mean | IQR | skewness | kurtosis | p00 | p01 | p05 | p10 | p20 | p25 | p30 | p40 | p50 | p60 | p70 | p75 | p80 | p90 | p95 | p99 | p100 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| feature_id | 6142140 | 0 | 1.510932e+07 | 6.036746e+06 | 2435.8078297 | 9.644721e+06 | -0.6648921 | -0.7973649 | 1810686.000 | 2299428.39000 | 3.668144e+06 | 4.736120e+06 | 8.471919e+06 | 1.024175e+07 | 1.278697e+07 | 1.538693e+07 | 1.700289e+07 | 1.821829e+07 | 1.933637e+07 | 1.988647e+07 | 2.068347e+07 | 2.191521e+07 | 2.234452e+07 | 2.274025e+07 | 22820741.000 |
| age | 6142140 | 0 | 1.343294e+02 | 8.024157e+01 | 0.0323772 | 8.800000e+01 | 0.9779711 | 1.1327092 | 0.000 | 10.00000 | 2.200000e+01 | 3.800000e+01 | 7.300000e+01 | 8.200000e+01 | 8.900000e+01 | 1.030000e+02 | 1.220000e+02 | 1.390000e+02 | 1.580000e+02 | 1.700000e+02 | 1.880000e+02 | 2.490000e+02 | 3.040000e+02 | 3.670000e+02 | 947.000 |
| vri_vol | 6142140 | 0 | 1.325707e+02 | 1.686811e+02 | 0.0680623 | 1.903580e+02 | 2.1607926 | 6.4878226 | 0.000 | 0.00000 | 0.000000e+00 | 0.000000e+00 | 3.328000e+00 | 8.192000e+00 | 1.554000e+01 | 3.773100e+01 | 7.081400e+01 | 1.133410e+02 | 1.664783e+02 | 1.985500e+02 | 2.355390e+02 | 3.430070e+02 | 4.626960e+02 | 7.803842e+02 | 2348.624 |
| l_vol | 6142140 | 0 | 1.278455e+02 | 1.692751e+02 | 0.0683020 | 1.838893e+02 | 2.2896282 | 7.3040882 | 0.000 | 0.00000 | 0.000000e+00 | 0.000000e+00 | 1.957380e+00 | 5.557100e+00 | 1.198660e+01 | 3.239386e+01 | 6.420535e+01 | 1.056077e+02 | 1.575780e+02 | 1.894464e+02 | 2.266294e+02 | 3.362262e+02 | 4.583625e+02 | 7.928209e+02 | 2700.065 |
| u_vol | 6142140 | 0 | 1.387145e+02 | 1.733656e+02 | 0.0699525 | 1.935548e+02 | 2.2149550 | 6.9264023 | 0.000 | 0.00000 | 0.000000e+00 | 0.000000e+00 | 5.746600e+00 | 1.180460e+01 | 1.989727e+01 | 4.390110e+01 | 7.836990e+01 | 1.206287e+02 | 1.734926e+02 | 2.053594e+02 | 2.419193e+02 | 3.502249e+02 | 4.747595e+02 | 8.149388e+02 | 2700.065 |
| vdyp_vol | 6142140 | 0 | 1.335848e+02 | 1.712624e+02 | 0.0691038 | 1.890707e+02 | 2.2462145 | 7.0755634 | 0.000 | 0.00000 | 0.000000e+00 | 0.000000e+00 | 4.090718e+00 | 8.858300e+00 | 1.612527e+01 | 3.852665e+01 | 7.163152e+01 | 1.135827e+02 | 1.661077e+02 | 1.979290e+02 | 2.347680e+02 | 3.436486e+02 | 4.666621e+02 | 8.031891e+02 | 2700.065 |
| est_age_spp1 | 6142139 | 1 | 1.215725e+02 | 7.784334e+01 | 0.0314095 | 9.000000e+01 | 0.9320963 | 1.0737184 | 1.000 | 2.00000 | 1.200000e+01 | 2.500000e+01 | 6.000000e+01 | 7.000000e+01 | 8.000000e+01 | 9.200000e+01 | 1.100000e+02 | 1.280000e+02 | 1.480000e+02 | 1.600000e+02 | 1.800000e+02 | 2.300000e+02 | 2.800000e+02 | 3.500000e+02 | 915.000 |
| vdyp_vol_est_age | 4768472 | 1373668 | 1.513132e+02 | 1.721776e+02 | 0.0788473 | 1.856472e+02 | 2.1940026 | 6.7933977 | 0.000 | 0.43316 | 2.833228e+00 | 6.435680e+00 | 1.925168e+01 | 2.814197e+01 | 3.870926e+01 | 6.456234e+01 | 9.637813e+01 | 1.350538e+02 | 1.840605e+02 | 2.137891e+02 | 2.489824e+02 | 3.603030e+02 | 4.913802e+02 | 8.196342e+02 | 2700.065 |
| diff_vdyp | 6142140 | 0 | 4.276689e+00 | 1.230657e+01 | 0.0049657 | 4.020357e+00 | 23.1306645 | 1405.4835331 | 0.000 | 0.00000 | 0.000000e+00 | 0.000000e+00 | 4.360000e-02 | 9.238880e-02 | 1.643014e-01 | 4.061800e-01 | 8.672440e-01 | 1.670316e+00 | 3.054901e+00 | 4.112746e+00 | 5.540973e+00 | 1.083846e+01 | 1.805625e+01 | 4.741201e+01 | 1786.505 |
| diff_vdyp_est_age | 4768472 | 1373668 | 2.111881e+01 | 3.979223e+01 | 0.0182225 | 2.218056e+01 | 5.7971005 | 63.5839205 | 0.000 | 0.01939 | 1.523700e-01 | 3.921000e-01 | 1.129720e+00 | 1.651600e+00 | 2.305090e+00 | 4.163980e+00 | 7.141825e+00 | 1.178884e+01 | 1.890183e+01 | 2.383216e+01 | 3.049928e+01 | 5.630750e+01 | 8.665497e+01 | 1.793640e+02 | 1844.900 |
| vol_dif | 6142140 | 0 | -1.014112e+00 | 1.298897e+01 | 0.0052410 | 1.779890e+00 | -16.0575109 | 1153.4680041 | -1786.505 | -40.69281 | -1.289980e+01 | -6.568179e+00 | -2.134842e+00 | -1.220650e+00 | -6.451400e-01 | -9.547600e-02 | 0.000000e+00 | 9.180000e-03 | 2.501020e-01 | 5.592400e-01 | 1.212572e+00 | 4.617402e+00 | 9.040567e+00 | 2.364802e+01 | 1771.781 |
## # A tibble: 1 × 10
## features mean_diff sd_diff se_mean diff_ge_5 pct_ge_5 diff_ge_10 pct_ge_10
## <int> <dbl> <dbl> <dbl> <int> <dbl> <int> <dbl>
## 1 6142140 -1.01 13.0 0.00524 1334337 21.7 677153 11.0
## # ℹ 2 more variables: diff_ge_20 <int>, pct_ge_20 <dbl>## # A tibble: 1 × 6
## features diff_ge_5 pct_ge_5 vdyp_gt_vri pct_vdyp_gt_vri vdyp_zero
## <int> <int> <dbl> <int> <dbl> <int>
## 1 597197 303969 50.9 379169 63.5 27680## # A tibble: 1 × 9
## features diff_ge_5 pct_ge_5 diff_ge_10 pct_ge_10 diff_ge_20 pct_ge_20
## <int> <int> <dbl> <int> <dbl> <int> <dbl>
## 1 5537525 1028709 18.6 458316 8.28 129640 2.34
## # ℹ 2 more variables: vri_zero_vdyp_positive <int>, vdyp_zero <int>Of the 6142140 features assessed the mean difference in volume between the VDYP dataset and the R1 is -1.01m3 (sd: 12.98), the standard error of the mean is 0.005m3 , 1334337 (21%) of the features have volume difference > +/- 5m3, 677153 (11%) features have volume difference > +/- 10m3 and 264097 (4%) features have volume difference > +/- 20m3.
The are a number of sources that drive divergence including approaches to interpolation (most common source for volume differences < 5m3), age cutoffs and layer/project table input error in both the R1 and vdyp_vol dataset. The provincial VDYP dataset has an age cutoff of 250 years effecting 3039699 (22%) of the 1334337 of the features have volume difference > +/- 5m3. 63% of the features > 250 years of age have vdyp volumes > the R1 volume suggesting the much older stands in the R1 are decrementing volume, driving the divergence. There are 27680 features that have ages > 250 and vdyp_vols == 0.
Of the features < 250 years of age 1028709 (18%) of the features have volume difference > +/- 5m3, 458316 (8%) features have volume difference > +/- 10m3 and 129640 (2%) features have volume difference > +/- 20m3. There 134276 features where the vri_vol == 0 and the vdyp_vol > 0 and 1099 features where the vdyp_vol == 0.
112262 features < 250 years have R1 yields more reflective of the layer table est_age_1 and 2 fields than the projected age (proj_age_1) suggesting an error has occured in the 2025 update.