Yield curves from Alberta Permanent Sample Plot (PSP) Data
Glen W Armstrong
Initialise
library(tidyverse)## -- Attaching packages --------------------------------------- tidyverse 1.3.2 --
## v ggplot2 3.4.0 v purrr 1.0.0
## v tibble 3.1.8 v dplyr 1.0.10
## v tidyr 1.2.1 v stringr 1.5.0
## v readr 2.1.3 v forcats 0.5.2
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()library(gt)Read in PSP data, and examine its structure. There seem to be problems with the opening number column, but we don’t need that so we will ignore.
#setwd("d:/gwa/AlbertaPSP_working/src")
pspCompiled <- read_csv('../data/compiled.csv')## Warning: One or more parsing issues, call `problems()` on your data frame for details,
## e.g.:
## dat <- vroom(...)
## problems(dat)## Rows: 193500 Columns: 119
## -- Column specification --------------------------------------------------------
## Delimiter: ","
## chr (21): company, fmu, fma, aspect, utm_zone, datum, natural_subregion, eco...
## dbl (91): company_plot_number, company_stand_number, establishment_year, est...
## lgl (7): sampling_unit_number, plot_comment, su, regen_decht, th_n_m, th_m,...
##
## i Use `spec()` to retrieve the full column specification for this data.
## i Specify the column types or set `show_col_types = FALSE` to quiet this message.problems(pspCompiled)## # A tibble: 160 x 5
## row col expected actual file
## <int> <int> <chr> <chr> <chr>
## 1 170622 13 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 2 170622 63 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 3 170623 13 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 4 170623 63 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 5 170624 13 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 6 170624 63 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 7 170625 13 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 8 170625 63 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 9 170626 13 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## 10 170626 63 a double 5210622265A C:/Users/gwa/Desktop/YieldCurveDevelopment~
## # ... with 150 more rowsstr(pspCompiled)## spc_tbl_ [193,500 x 119] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ company : chr [1:193500] "GOA" "GOA" "GOA" "GOA" ...
## $ company_plot_number : num [1:193500] 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 ...
## $ company_stand_number: num [1:193500] 1 1 1 1 1 1 1 1 1 1 ...
## $ establishment_year : num [1:193500] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960 ...
## $ establishment_month : num [1:193500] 7 7 7 7 7 7 7 7 7 7 ...
## $ establishment_day : num [1:193500] 23 23 23 23 23 23 23 23 23 23 ...
## $ fmu : chr [1:193500] "W14" "W14" "W14" "W14" ...
## $ fma : chr [1:193500] "FMA7500020" "FMA7500020" "FMA7500020" "FMA7500020" ...
## $ ats_township : num [1:193500] 62 62 62 62 62 62 62 62 62 62 ...
## $ ats_range : num [1:193500] 19 19 19 19 19 19 19 19 19 19 ...
## $ ats_meridian : num [1:193500] 5 5 5 5 5 5 5 5 5 5 ...
## $ ats_section : num [1:193500] 31 31 31 31 31 31 31 31 31 31 ...
## $ opening_number : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ sampling_unit_number: logi [1:193500] NA NA NA NA NA NA ...
## $ topographic_position: num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ elevation : num [1:193500] 811 811 811 811 811 811 811 811 811 811 ...
## $ slope : num [1:193500] 2 2 2 2 2 2 2 2 2 2 ...
## $ aspect : chr [1:193500] "SW" "SW" "SW" "SW" ...
## $ x_coord : num [1:193500] 510771 510771 510771 510771 510771 ...
## $ y_coord : num [1:193500] 6028487 6028487 6028487 6028487 6028487 ...
## $ utm_zone : chr [1:193500] "UTM11" "UTM11" "UTM11" "UTM11" ...
## $ datum : chr [1:193500] "NAD83" "NAD83" "NAD83" "NAD83" ...
## $ latitude : num [1:193500] 54.4 54.4 54.4 54.4 54.4 ...
## $ longitude : num [1:193500] -117 -117 -117 -117 -117 ...
## $ natural_subregion : chr [1:193500] "CM" "CM" "CM" "CM" ...
## $ ecosite_guide : chr [1:193500] "N" "N" "N" "N" ...
## $ ecosite : chr [1:193500] NA NA NA NA ...
## $ ecosite_phase : chr [1:193500] NA NA NA NA ...
## $ plot_comment : logi [1:193500] NA NA NA NA NA NA ...
## $ shared : num [1:193500] 1 1 1 1 1 1 1 1 1 1 ...
## $ strata : chr [1:193500] "HwSw" "HwSw" "HwSw" "HwSw" ...
## $ substrata : chr [1:193500] "Pb" "Pb" "Pb" "Pb" ...
## $ est_yr : num [1:193500] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960 ...
## $ mmt_num : num [1:193500] 1 1 1 1 1 1 1 1 1 1 ...
## $ mmt_yr : num [1:193500] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960 ...
## $ mmt_mo : num [1:193500] 7 7 7 7 7 7 7 7 7 7 ...
## $ age : num [1:193500] 81 81 81 81 81 81 81 81 81 81 ...
## $ stand_type : num [1:193500] 1 1 1 1 1 1 1 1 1 1 ...
## $ scale : chr [1:193500] "con_dec" "con_dec" "con_dec" "species" ...
## $ species : chr [1:193500] "CO" "DE" "TO" "AW" ...
## $ sphRegenOnly : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ sphBH : num [1:193500] 691.7 306.3 998 108.7 39.5 ...
## $ sphD15 : num [1:193500] 691.7 306.3 998 108.7 39.5 ...
## $ sphD91 : num [1:193500] 652.2 306.3 958.5 108.7 39.5 ...
## $ ba : num [1:193500] 20.814 22.01 42.825 8.801 0.625 ...
## $ vol_0000 : num [1:193500] 167.04 190.07 357.12 86.6 3.55 ...
## $ vol_1307 : num [1:193500] 153.01 181.23 334.24 83.33 2.74 ...
## $ vol_1510 : num [1:193500] 141.99 177.84 319.83 82.32 2.02 ...
## $ biomass : num [1:193500] 79398 89934 169331 43342 2819 ...
## $ carbon : num [1:193500] 39699 44967 84666 21671 1409 ...
## $ topht_n : num [1:193500] NA NA NA 10 NA NA 8 NA NA 10 ...
## $ topht : num [1:193500] NA NA NA 23.9 13.7 ...
## $ topht_stat : chr [1:193500] NA NA NA "MP" ...
## $ age_n : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ age_tot : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ age_bh : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ age_max : num [1:193500] 81.4 81.4 81.4 81.4 81.4 ...
## $ age_harv : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ age_fire : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ age_avi : num [1:193500] 90 90 90 90 90 90 90 90 90 90 ...
## $ si_n : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ si_bh : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ opening : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ su : logi [1:193500] NA NA NA NA NA NA ...
## $ topo : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ elev : num [1:193500] 811 811 811 811 811 811 811 811 811 811 ...
## $ utm : chr [1:193500] "UTM11" "UTM11" "UTM11" "UTM11" ...
## $ natsub : chr [1:193500] "CM" "CM" "CM" "CM" ...
## $ ecoguide : chr [1:193500] "N" "N" "N" "N" ...
## $ eco : chr [1:193500] NA NA NA NA ...
## $ ecophs : chr [1:193500] NA NA NA NA ...
## $ stand_origin : chr [1:193500] "N" "N" "N" "N" ...
## $ plot_type : num [1:193500] 3 3 3 3 3 3 3 3 3 3 ...
## $ tree_area : num [1:193500] 1012 1012 1012 1012 1012 ...
## $ tree_tag : num [1:193500] 9.1 9.1 9.1 9.1 9.1 9.1 9.1 9.1 9.1 9.1 ...
## $ sap_area : num [1:193500] 1012 1012 1012 1012 1012 ...
## $ sap_dbh : num [1:193500] 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 ...
## $ sap_ht : num [1:193500] 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 ...
## $ regen_area : num [1:193500] 0 0 0 0 0 0 0 0 0 0 ...
## $ regen_conht : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ regen_decht : logi [1:193500] NA NA NA NA NA NA ...
## $ regen_n : num [1:193500] 0 0 0 0 0 0 0 0 0 0 ...
## $ regendone : num [1:193500] 0 0 0 0 0 0 0 0 0 0 ...
## $ fire_orig : num [1:193500] NA NA NA NA NA NA NA NA NA NA ...
## $ target : num [1:193500] NA NA NA 10 10 NA 10 NA NA 10 ...
## $ _FREQ_ : num [1:193500] NA NA NA 10 4 NA 6 NA NA 10 ...
## $ th_aa : num [1:193500] NA NA NA 23.9 13.7 ...
## $ th_n_a : num [1:193500] NA NA NA 10 4 NA 8 NA NA 10 ...
## $ th_mp : num [1:193500] NA NA NA 23.9 NA ...
## $ th_n_mp : num [1:193500] NA NA NA 10 NA NA 8 NA NA 10 ...
## $ th_mm : num [1:193500] NA NA NA 21.7 NA ...
## $ th_n_mm : num [1:193500] NA NA NA 3 NA NA NA NA NA 2 ...
## $ th_ms : num [1:193500] NA NA NA 21.7 NA ...
## $ th_n_ms : num [1:193500] NA NA NA 3 0 NA 0 NA NA 1 ...
## $ th_mpe : num [1:193500] NA NA NA 23.9 13.7 ...
## $ th_n_mpe : num [1:193500] NA NA NA 10 4 NA 6 NA NA 10 ...
## $ th_n_m : logi [1:193500] NA NA NA NA NA NA ...
## $ th_m : logi [1:193500] NA NA NA NA NA NA ...
## $ th_n_aa : logi [1:193500] NA NA NA NA NA NA ...
## [list output truncated]
## - attr(*, "spec")=
## .. cols(
## .. company = col_character(),
## .. company_plot_number = col_double(),
## .. company_stand_number = col_double(),
## .. establishment_year = col_double(),
## .. establishment_month = col_double(),
## .. establishment_day = col_double(),
## .. fmu = col_character(),
## .. fma = col_character(),
## .. ats_township = col_double(),
## .. ats_range = col_double(),
## .. ats_meridian = col_double(),
## .. ats_section = col_double(),
## .. opening_number = col_double(),
## .. sampling_unit_number = col_logical(),
## .. topographic_position = col_double(),
## .. elevation = col_double(),
## .. slope = col_double(),
## .. aspect = col_character(),
## .. x_coord = col_double(),
## .. y_coord = col_double(),
## .. utm_zone = col_character(),
## .. datum = col_character(),
## .. latitude = col_double(),
## .. longitude = col_double(),
## .. natural_subregion = col_character(),
## .. ecosite_guide = col_character(),
## .. ecosite = col_character(),
## .. ecosite_phase = col_character(),
## .. plot_comment = col_logical(),
## .. shared = col_double(),
## .. strata = col_character(),
## .. substrata = col_character(),
## .. est_yr = col_double(),
## .. mmt_num = col_double(),
## .. mmt_yr = col_double(),
## .. mmt_mo = col_double(),
## .. age = col_double(),
## .. stand_type = col_double(),
## .. scale = col_character(),
## .. species = col_character(),
## .. sphRegenOnly = col_double(),
## .. sphBH = col_double(),
## .. sphD15 = col_double(),
## .. sphD91 = col_double(),
## .. ba = col_double(),
## .. vol_0000 = col_double(),
## .. vol_1307 = col_double(),
## .. vol_1510 = col_double(),
## .. biomass = col_double(),
## .. carbon = col_double(),
## .. topht_n = col_double(),
## .. topht = col_double(),
## .. topht_stat = col_character(),
## .. age_n = col_double(),
## .. age_tot = col_double(),
## .. age_bh = col_double(),
## .. age_max = col_double(),
## .. age_harv = col_double(),
## .. age_fire = col_double(),
## .. age_avi = col_double(),
## .. si_n = col_double(),
## .. si_bh = col_double(),
## .. opening = col_double(),
## .. su = col_logical(),
## .. topo = col_double(),
## .. elev = col_double(),
## .. utm = col_character(),
## .. natsub = col_character(),
## .. ecoguide = col_character(),
## .. eco = col_character(),
## .. ecophs = col_character(),
## .. stand_origin = col_character(),
## .. plot_type = col_double(),
## .. tree_area = col_double(),
## .. tree_tag = col_double(),
## .. sap_area = col_double(),
## .. sap_dbh = col_double(),
## .. sap_ht = col_double(),
## .. regen_area = col_double(),
## .. regen_conht = col_double(),
## .. regen_decht = col_logical(),
## .. regen_n = col_double(),
## .. regendone = col_double(),
## .. fire_orig = col_double(),
## .. target = col_double(),
## .. `_FREQ_` = col_double(),
## .. th_aa = col_double(),
## .. th_n_a = col_double(),
## .. th_mp = col_double(),
## .. th_n_mp = col_double(),
## .. th_mm = col_double(),
## .. th_n_mm = col_double(),
## .. th_ms = col_double(),
## .. th_n_ms = col_double(),
## .. th_mpe = col_double(),
## .. th_n_mpe = col_double(),
## .. th_n_m = col_logical(),
## .. th_m = col_logical(),
## .. th_n_aa = col_logical(),
## .. totba = col_double(),
## .. totregden = col_double(),
## .. conba = col_double(),
## .. conregden = col_double(),
## .. decba = col_double(),
## .. decregden = col_double(),
## .. pinba = col_double(),
## .. pinregden = col_double(),
## .. sfba = col_double(),
## .. sfregden = col_double(),
## .. sbba = col_double(),
## .. sbregden = col_double(),
## .. pjba = col_double(),
## .. pjregden = col_double(),
## .. firba = col_double(),
## .. firregden = col_double(),
## .. bwba = col_double(),
## .. bwregden = col_double(),
## .. pbba = col_double(),
## .. pbregden = col_double()
## .. )
## - attr(*, "problems")=<externalptr>There’s a whole bunch of stuff we don’t need. Keep only relevant info.
psp <- pspCompiled |> select(company_plot_number, mmt_num, mmt_yr, species,
vol_1510, elevation, latitude, longitude, natsub,
strata, scale,age)
str(psp)## tibble [193,500 x 12] (S3: tbl_df/tbl/data.frame)
## $ company_plot_number: num [1:193500] 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 1.1 ...
## $ mmt_num : num [1:193500] 1 1 1 1 1 1 1 1 1 1 ...
## $ mmt_yr : num [1:193500] 1960 1960 1960 1960 1960 1960 1960 1960 1960 1960 ...
## $ species : chr [1:193500] "CO" "DE" "TO" "AW" ...
## $ vol_1510 : num [1:193500] 141.99 177.84 319.83 82.32 2.02 ...
## $ elevation : num [1:193500] 811 811 811 811 811 811 811 811 811 811 ...
## $ latitude : num [1:193500] 54.4 54.4 54.4 54.4 54.4 ...
## $ longitude : num [1:193500] -117 -117 -117 -117 -117 ...
## $ natsub : chr [1:193500] "CM" "CM" "CM" "CM" ...
## $ strata : chr [1:193500] "HwSw" "HwSw" "HwSw" "HwSw" ...
## $ scale : chr [1:193500] "con_dec" "con_dec" "con_dec" "species" ...
## $ age : num [1:193500] 81 81 81 81 81 81 81 81 81 81 ...Filtering and mutating.
psp <- psp |>
filter(scale =="con_dec") |>
filter(natsub == "CM" | natsub == "DM")
dim(pspCompiled)## [1] 193500 119dim(psp)## [1] 5517 12head(psp)## # A tibble: 6 x 12
## company~1 mmt_num mmt_yr species vol_1~2 eleva~3 latit~4 longi~5 natsub strata
## <dbl> <dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <chr> <chr>
## 1 1.1 1 1960 CO 142. 811 54.4 -117. CM HwSw
## 2 1.1 1 1960 DE 178. 811 54.4 -117. CM HwSw
## 3 1.1 1 1960 TO 320. 811 54.4 -117. CM HwSw
## 4 1.1 2 1965 CO 160. 811 54.4 -117. CM HwSw
## 5 1.1 2 1965 DE 187. 811 54.4 -117. CM HwSw
## 6 1.1 2 1965 TO 348. 811 54.4 -117. CM HwSw
## # ... with 2 more variables: scale <chr>, age <dbl>, and abbreviated variable
## # names 1: company_plot_number, 2: vol_1510, 3: elevation, 4: latitude,
## # 5: longitudetail(psp)## # A tibble: 6 x 12
## company~1 mmt_num mmt_yr species vol_1~2 eleva~3 latit~4 longi~5 natsub strata
## <dbl> <dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <chr> <chr>
## 1 938. 1 1990 CO 70.2 587 55.2 -114. CM Sb
## 2 938. 1 1990 DE 4.22 587 55.2 -114. CM Sb
## 3 938. 1 1990 TO 74.4 587 55.2 -114. CM Sb
## 4 938. 2 1995 CO 99.6 587 55.2 -114. CM Sb
## 5 938. 2 1995 DE 5.33 587 55.2 -114. CM Sb
## 6 938. 2 1995 TO 105. 587 55.2 -114. CM Sb
## # ... with 2 more variables: scale <chr>, age <dbl>, and abbreviated variable
## # names 1: company_plot_number, 2: vol_1510, 3: elevation, 4: latitude,
## # 5: longitudeThe plots are broken down by subplot. We don’t want that. Aggregate!
psp <- psp |>
mutate(plotnum = floor(company_plot_number)) |>
mutate(subplotnum = (company_plot_number - plotnum) * 10)
head(psp)## # A tibble: 6 x 14
## company~1 mmt_num mmt_yr species vol_1~2 eleva~3 latit~4 longi~5 natsub strata
## <dbl> <dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <chr> <chr>
## 1 1.1 1 1960 CO 142. 811 54.4 -117. CM HwSw
## 2 1.1 1 1960 DE 178. 811 54.4 -117. CM HwSw
## 3 1.1 1 1960 TO 320. 811 54.4 -117. CM HwSw
## 4 1.1 2 1965 CO 160. 811 54.4 -117. CM HwSw
## 5 1.1 2 1965 DE 187. 811 54.4 -117. CM HwSw
## 6 1.1 2 1965 TO 348. 811 54.4 -117. CM HwSw
## # ... with 4 more variables: scale <chr>, age <dbl>, plotnum <dbl>,
## # subplotnum <dbl>, and abbreviated variable names 1: company_plot_number,
## # 2: vol_1510, 3: elevation, 4: latitude, 5: longitudetail(psp)## # A tibble: 6 x 14
## company~1 mmt_num mmt_yr species vol_1~2 eleva~3 latit~4 longi~5 natsub strata
## <dbl> <dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <chr> <chr>
## 1 938. 1 1990 CO 70.2 587 55.2 -114. CM Sb
## 2 938. 1 1990 DE 4.22 587 55.2 -114. CM Sb
## 3 938. 1 1990 TO 74.4 587 55.2 -114. CM Sb
## 4 938. 2 1995 CO 99.6 587 55.2 -114. CM Sb
## 5 938. 2 1995 DE 5.33 587 55.2 -114. CM Sb
## 6 938. 2 1995 TO 105. 587 55.2 -114. CM Sb
## # ... with 4 more variables: scale <chr>, age <dbl>, plotnum <dbl>,
## # subplotnum <dbl>, and abbreviated variable names 1: company_plot_number,
## # 2: vol_1510, 3: elevation, 4: latitude, 5: longitudedim(psp)## [1] 5517 14summary(psp)## company_plot_number mmt_num mmt_yr species
## Min. : 1.1 Min. :1.000 Min. :1960 Length:5517
## 1st Qu.: 256.2 1st Qu.:1.000 1st Qu.:1980 Class :character
## Median : 337.4 Median :3.000 Median :1990 Mode :character
## Mean :2088.2 Mean :2.993 Mean :1989
## 3rd Qu.:3013.1 3rd Qu.:4.000 3rd Qu.:1998
## Max. :9033.1 Max. :9.000 Max. :2019
##
## vol_1510 elevation latitude longitude
## Min. : 0.00 Min. :302 Min. :53.03 Min. :-119.9
## 1st Qu.: 0.00 1st Qu.:381 1st Qu.:54.78 1st Qu.:-117.6
## Median : 89.44 Median :624 Median :55.33 Median :-116.8
## Mean :135.63 Mean :597 Mean :56.08 Mean :-115.9
## 3rd Qu.:245.47 3rd Qu.:754 3rd Qu.:57.82 3rd Qu.:-114.3
## Max. :786.14 Max. :898 Max. :59.19 Max. :-110.1
## NA's :192
## natsub strata scale age
## Length:5517 Length:5517 Length:5517 Min. :-16.00
## Class :character Class :character Class :character 1st Qu.: 27.00
## Mode :character Mode :character Mode :character Median : 82.00
## Mean : 79.49
## 3rd Qu.:123.00
## Max. :246.00
## NA's :663
## plotnum subplotnum
## Min. : 1 Min. :1.000
## 1st Qu.: 256 1st Qu.:1.000
## Median : 337 Median :1.000
## Mean :2088 Mean :1.794
## 3rd Qu.:3013 3rd Qu.:3.000
## Max. :9033 Max. :4.000
## Let’s aggregate to the plot level.
psp_agg <- psp |>
arrange(plotnum,mmt_num,species) |>
group_by(plotnum,mmt_num,species) |>
summarise(latitude = mean(latitude), longitude = mean(longitude),
elevation = mean(elevation), vol_1510 = mean(vol_1510),
age = mean(age), mmt_yr = mean(mmt_yr))## `summarise()` has grouped output by 'plotnum', 'mmt_num'. You can override
## using the `.groups` argument.head(psp_agg)## # A tibble: 6 x 9
## # Groups: plotnum, mmt_num [2]
## plotnum mmt_num species latitude longitude elevation vol_1510 age mmt_yr
## <dbl> <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 CO 54.4 -117. 811. 184. 80 1960
## 2 1 1 DE 54.4 -117. 811. 128. 80 1960
## 3 1 1 TO 54.4 -117. 811. 312. 80 1960
## 4 1 2 CO 54.4 -117. 811. 202. 85 1965
## 5 1 2 DE 54.4 -117. 811. 132. 85 1965
## 6 1 2 TO 54.4 -117. 811. 333. 85 1965dim(psp_agg)## [1] 3333 9Pivot so that volumes by species group are all in one row.
psp_pivot <- psp_agg |>
pivot_wider(names_from = species, values_from= vol_1510)
dim(psp_pivot)## [1] 1111 10head(psp_pivot)## # A tibble: 6 x 10
## # Groups: plotnum, mmt_num [6]
## plotnum mmt_num latitude longitude elevation age mmt_yr CO DE TO
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 54.4 -117. 811. 80 1960 184. 128. 312.
## 2 1 2 54.4 -117. 811. 85 1965 202. 132. 333.
## 3 1 3 54.4 -117. 811. 88 1968 216. 137. 353.
## 4 1 4 54.4 -117. 811. 100 1980 231. 131. 361.
## 5 1 5 54.4 -117. 811. 110 1990 224. 145. 369.
## 6 1 6 54.4 -117. 811. 122 2002 177. 92.1 269.summary(psp_pivot)## plotnum mmt_num latitude longitude
## Min. : 1 Min. :1.000 Min. :53.03 Min. :-119.9
## 1st Qu.: 376 1st Qu.:1.000 1st Qu.:55.07 1st Qu.:-117.4
## Median : 641 Median :3.000 Median :55.64 Median :-115.3
## Mean :3309 Mean :3.167 Mean :56.13 Mean :-115.0
## 3rd Qu.:8185 3rd Qu.:4.000 3rd Qu.:57.19 3rd Qu.:-113.1
## Max. :9033 Max. :9.000 Max. :59.19 Max. :-110.1
##
## elevation age mmt_yr CO
## Min. :302.2 Min. :-16.00 Min. :1960 Min. : 0.0000
## 1st Qu.:400.0 1st Qu.: 11.50 1st Qu.:1987 1st Qu.: 0.0000
## Median :608.0 Median : 34.00 Median :1993 Median : 0.9925
## Mean :567.7 Mean : 52.12 Mean :1993 Mean : 67.6428
## 3rd Qu.:644.0 3rd Qu.: 84.12 3rd Qu.:2001 3rd Qu.:117.4875
## Max. :898.0 Max. :211.50 Max. :2019 Max. :488.0192
## NA's :49 NA's :192
## DE TO
## Min. : 0.000 Min. : 0.00
## 1st Qu.: 0.000 1st Qu.: 0.00
## Median : 5.751 Median : 55.67
## Mean : 60.714 Mean :128.36
## 3rd Qu.: 92.172 3rd Qu.:247.85
## Max. :489.629 Max. :595.12
## Let’s recover the stratum information. It’s possible that subplots within a plot have different strata. Use the first one encountered.
psp_strat <- psp |>
arrange(plotnum,mmt_num) |>
group_by(plotnum,mmt_num) |>
summarise(stratum = first(strata))## `summarise()` has grouped output by 'plotnum'. You can override using the
## `.groups` argument.head(psp_strat)## # A tibble: 6 x 3
## # Groups: plotnum [1]
## plotnum mmt_num stratum
## <dbl> <dbl> <chr>
## 1 1 1 HwSw
## 2 1 2 HwSw
## 3 1 3 HwSw
## 4 1 4 HwSw
## 5 1 5 HwSw
## 6 1 6 HwSwdim(psp_strat)## [1] 1111 3summary(psp_strat)## plotnum mmt_num stratum
## Min. : 1 Min. :1.000 Length:1111
## 1st Qu.: 376 1st Qu.:1.000 Class :character
## Median : 641 Median :3.000 Mode :character
## Mean :3309 Mean :3.167
## 3rd Qu.:8185 3rd Qu.:4.000
## Max. :9033 Max. :9.000Add stratum back to aggregated data
mergeCols <- c("plotnum","mmt_num")
psp_merged <- merge(psp_pivot, psp_strat, by = mergeCols)
dim(psp_merged)## [1] 1111 11head(psp_merged)## plotnum mmt_num latitude longitude elevation age mmt_yr CO DE
## 1 1 1 54.40416 -116.8351 811.25 80 1960 184.1957 127.69273
## 2 1 2 54.40416 -116.8351 811.25 85 1965 201.5248 131.80729
## 3 1 3 54.40416 -116.8351 811.25 88 1968 215.6712 136.90179
## 4 1 4 54.40416 -116.8351 811.25 100 1980 230.7656 130.50700
## 5 1 5 54.40416 -116.8351 811.25 110 1990 223.7412 144.97638
## 6 1 6 54.40416 -116.8351 811.25 122 2002 176.8404 92.10343
## TO stratum
## 1 311.8884 HwSw
## 2 333.3321 HwSw
## 3 352.5730 HwSw
## 4 361.2726 HwSw
## 5 368.7176 HwSw
## 6 268.9439 HwSwsummary(psp_merged)## plotnum mmt_num latitude longitude
## Min. : 1 Min. :1.000 Min. :53.03 Min. :-119.9
## 1st Qu.: 376 1st Qu.:1.000 1st Qu.:55.07 1st Qu.:-117.4
## Median : 641 Median :3.000 Median :55.64 Median :-115.3
## Mean :3309 Mean :3.167 Mean :56.13 Mean :-115.0
## 3rd Qu.:8185 3rd Qu.:4.000 3rd Qu.:57.19 3rd Qu.:-113.1
## Max. :9033 Max. :9.000 Max. :59.19 Max. :-110.1
##
## elevation age mmt_yr CO
## Min. :302.2 Min. :-16.00 Min. :1960 Min. : 0.0000
## 1st Qu.:400.0 1st Qu.: 11.50 1st Qu.:1987 1st Qu.: 0.0000
## Median :608.0 Median : 34.00 Median :1993 Median : 0.9925
## Mean :567.7 Mean : 52.12 Mean :1993 Mean : 67.6428
## 3rd Qu.:644.0 3rd Qu.: 84.12 3rd Qu.:2001 3rd Qu.:117.4875
## Max. :898.0 Max. :211.50 Max. :2019 Max. :488.0192
## NA's :49 NA's :192
## DE TO stratum
## Min. : 0.000 Min. : 0.00 Length:1111
## 1st Qu.: 0.000 1st Qu.: 0.00 Class :character
## Median : 5.751 Median : 55.67 Mode :character
## Mean : 60.714 Mean :128.36
## 3rd Qu.: 92.172 3rd Qu.:247.85
## Max. :489.629 Max. :595.12
## Get rid of records with negative ages and rename volume columns.
psp_merged <- psp_merged |>
filter(age >= 0) |>
rename(vol_tot = TO, vol_con = CO, vol_dec = DE)
summary(psp_merged)## plotnum mmt_num latitude longitude
## Min. : 1.0 Min. :1.000 Min. :53.03 Min. :-119.4
## 1st Qu.: 349.0 1st Qu.:2.000 1st Qu.:55.08 1st Qu.:-117.4
## Median : 681.5 Median :3.000 Median :55.64 Median :-115.3
## Mean :3629.7 Mean :3.412 Mean :56.21 Mean :-114.9
## 3rd Qu.:9006.0 3rd Qu.:5.000 3rd Qu.:57.82 3rd Qu.:-113.0
## Max. :9033.0 Max. :9.000 Max. :59.19 Max. :-110.1
##
## elevation age mmt_yr vol_con
## Min. :303 Min. : 0.00 Min. :1960 Min. : 0.0
## 1st Qu.:381 1st Qu.: 12.00 1st Qu.:1987 1st Qu.: 0.0
## Median :608 Median : 35.00 Median :1992 Median : 0.0
## Mean :562 Mean : 53.26 Mean :1992 Mean : 64.7
## 3rd Qu.:640 3rd Qu.: 85.00 3rd Qu.:2001 3rd Qu.:102.8
## Max. :898 Max. :211.50 Max. :2019 Max. :488.0
## NA's :10
## vol_dec vol_tot stratum
## Min. : 0.000 Min. : 0.00 Length:902
## 1st Qu.: 0.000 1st Qu.: 0.00 Class :character
## Median : 1.803 Median : 16.83 Mode :character
## Mean : 50.393 Mean :115.09
## 3rd Qu.: 76.335 3rd Qu.:234.82
## Max. :489.629 Max. :595.12
## Analysis
Let’s find number of unique combinations of plot and stratum.
psp_merged |>
select(plotnum,stratum) |>
distinct(plotnum,stratum) |>
group_by(stratum) |>
summarise(count=n())## # A tibble: 7 x 2
## stratum count
## <chr> <int>
## 1 Hw 68
## 2 HwSw 41
## 3 Pl 21
## 4 PlHw 8
## 5 Sb 4
## 6 Sw 23
## 7 SwHw 30Decision We will use only two yield strata, mixed (Sw, SwHw, HwSw, and Hw), and pine-black spruce (Pl, PlHw, and Sb).
Graphical summary of PSP data
ggplot(psp_merged, aes(x = age, y = vol_tot, colour = stratum, group = plotnum)) +
geom_line() +
ggtitle("Total Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
ggplot(psp_merged, aes(x = age, y = vol_con, colour = stratum, group = plotnum)) +
geom_line() +
ggtitle("Softwood Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
ggplot(psp_merged, aes(x = age, y = vol_dec, colour = stratum, group = plotnum)) +
geom_line() +
ggtitle("Hardwood Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
Replot aggregated strata
psp_merged <- psp_merged |>
mutate(stratum2 = recode(stratum,
"Hw" = "Mixed",
"HwSw" = "Mixed",
"Sw" = "Mixed",
"SwHw" = "Mixed",
"Pl" = "PlSb",
"PlHw" = "PlSb",
"Sb" = "PlSb")) |>
arrange(stratum2,age)Save psp_merged so that it can be read in later.
save(psp_merged, file = "psp_merged.RData")ggplot(psp_merged, aes(x = age, y = vol_tot, colour = stratum2, group = plotnum)) +
geom_line() +
ggtitle("Total Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
ggplot(psp_merged, aes(x = age, y = vol_con, colour = stratum2, group = plotnum)) +
geom_line() +
ggtitle("Softwood Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
ggplot(psp_merged, aes(x = age, y = vol_dec, colour = stratum2, group = plotnum)) +
geom_line() +
ggtitle("Hardwood Volume") +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)")
Estimation (non-linear regression)
We will estimate total and softwood yields using Chapman-Richards equation, and calculate hardwood yields as the difference between total and softwood yields. The Chapman-Richards equation is below, where \(t\) represents stand age in years.
\[ v(t) = b_0 (1 - exp(-b_1 t))^{b_2} \] ### Mixed stratum
mixed <- psp_merged |>
filter(stratum2 == "Mixed" )
x <- mixed$age
y <- mixed$vol_tot
m <- nls(y ~ b0 * (1 - exp(-b1 * x))^b2, start = list(b0=300,b1=0.027,b2=4.0))
summary(m)##
## Formula: y ~ b0 * (1 - exp(-b1 * x))^b2
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## b0 3.464e+02 4.880e+00 70.982 < 2e-16 ***
## b1 5.558e-02 4.029e-03 13.796 < 2e-16 ***
## b2 1.730e+01 3.713e+00 4.661 3.7e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 54.65 on 789 degrees of freedom
##
## Number of iterations to convergence: 12
## Achieved convergence tolerance: 3.104e-06plot(x,y,type='p',ylim=c(0,500))
par(new=TRUE)
plot(x,predict(m),col="red",type = "l", axes = FALSE, ylim = c(0,500))
coef(m)## b0 b1 b2
## 346.39508943 0.05558123 17.30446221coeff_table <- tibble_row(stratum="mixed", type="tot",b0=coef(m)[1],b1=coef(m)[2],b2=coef(m)[3])mixed <- psp_merged |>
filter(stratum2 == "Mixed")
y <- mixed$vol_con
x <- mixed$age
m <- nls(y~b0*(1-exp(-b1*x))^b2,start=list(b0=300,b1=0.027,b2=4.0))
plot(x,y,type='p',ylim=c(0,500))
par(new=TRUE)
plot(x,predict(m),col="red",type='l',axes=FALSE,ylim=c(0,500))
summary(m)##
## Formula: y ~ b0 * (1 - exp(-b1 * x))^b2
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## b0 2.801e+02 9.250e+00 30.284 <2e-16 ***
## b1 3.984e-02 4.145e-03 9.611 <2e-16 ***
## b2 2.535e+01 8.323e+00 3.046 0.0024 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 51.35 on 789 degrees of freedom
##
## Number of iterations to convergence: 10
## Achieved convergence tolerance: 2.441e-06coeff_table <- coeff_table |> add_row(stratum="mixed", type="con",b0=coef(m)[1],b1=coef(m)[2],b2=coef(m)[3])Pine/black spruce stratum
Total volume
plsb <- psp_merged |>
filter(stratum2 == "PlSb")
y <- plsb$vol_tot
x <- plsb$age
b0 <- 300
m <- nls(y~b0*(1-exp(-b1*x))^b2,start=list(b1=0.027,b2=4.0),control=nls.control(maxiter = 1000))
#m <- nls(y~b0*(1-exp(-b1*x))^b2,start=list(b0=450,b1=0.017702,b2=4.555790))
plot(x,y,type='p',ylim=c(0,500))
par(new=TRUE)
plot(x,predict(m),col="red",type='l',axes=FALSE,ylim=c(0,500))
summary(m)##
## Formula: y ~ b0 * (1 - exp(-b1 * x))^b2
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## b1 0.02755 0.00446 6.178 1.17e-08 ***
## b2 7.06801 2.39056 2.957 0.00382 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 64.2 on 108 degrees of freedom
##
## Number of iterations to convergence: 9
## Achieved convergence tolerance: 3.34e-06coeff_table <- coeff_table |> add_row(stratum="plsb", type="tot",b0,b1=coef(m)[1],b2=coef(m)[2])Conifer volume
Pine/black spruce stratum
plsb <- psp_merged |>
filter(stratum2 == "PlSb")
y <- plsb$vol_con
x <- plsb$age
b0 <- 300
m <- nls(y~b0*(1-exp(-b1*x))^b2,start=list(b1=0.027,b2=4.0),
control=nls.control(maxiter = 1000))
#m <- nls(y~b0*(1-exp(-b1*x))^b2,start=list(b0=450,b1=0.017702,b2=4.555790))
plot(x,y,type='p',ylim=c(0,500))
par(new=TRUE)
plot(x,predict(m),col="red",type='l',axes=FALSE,ylim=c(0,500))
summary(m)##
## Formula: y ~ b0 * (1 - exp(-b1 * x))^b2
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## b1 0.019011 0.002879 6.604 1.55e-09 ***
## b2 4.371952 1.073640 4.072 8.91e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 48.89 on 108 degrees of freedom
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 8.662e-06coeff_table <- coeff_table |> add_row(stratum="plsb", type="con",b0,b1=coef(m)[1],b2=coef(m)[2])Print out coefficients table
coeff_table## # A tibble: 4 x 5
## stratum type b0 b1 b2
## <chr> <chr> <dbl> <dbl> <dbl>
## 1 mixed tot 346. 0.0556 17.3
## 2 mixed con 280. 0.0398 25.4
## 3 plsb tot 300 0.0276 7.07
## 4 plsb con 300 0.0190 4.37Create a function to calculate volume (Chapman-Richards)
cr_volume <- function(b0,b1,b2,age){
vol <- b0 * (1 - exp(-b1 * age))^b2
return(vol)
}
cr_volume(364.4,0.05558,17.3,100)## [1] 340.8376Create a table of yields for mixedwood plots
age <- 1:40 * 5
mixed_table <- tibble(age, totvol = cr_volume(364.4,0.05558,17.3,age),
convol = cr_volume(280.1,0.03984,25.35,age),
decvol = totvol - convol)
mixed_table## # A tibble: 40 x 4
## age totvol convol decvol
## <dbl> <dbl> <dbl> <dbl>
## 1 5 8.34e-9 4.04e-17 8.34e-9
## 2 10 1.44e-4 1.57e-10 1.44e-4
## 3 15 1.90e-2 4.50e- 7 1.90e-2
## 4 20 3.66e-1 7.08e- 5 3.66e-1
## 5 25 2.56e+0 2.35e- 3 2.56e+0
## 6 30 9.77e+0 3.01e- 2 9.74e+0
## 7 35 2.53e+1 2.04e- 1 2.51e+1
## 8 40 5.02e+1 8.84e- 1 4.93e+1
## 9 45 8.29e+1 2.77e+ 0 8.02e+1
## 10 50 1.20e+2 6.80e+ 0 1.13e+2
## # ... with 30 more rowsmixed_table2 <-
mixed_table |>
pivot_longer(cols=c("totvol","convol","decvol"), names_to = "type",
values_to = "volume") |>
arrange(type,age)
mixed_table2## # A tibble: 120 x 3
## age type volume
## <dbl> <chr> <dbl>
## 1 5 convol 4.04e-17
## 2 10 convol 1.57e-10
## 3 15 convol 4.50e- 7
## 4 20 convol 7.08e- 5
## 5 25 convol 2.35e- 3
## 6 30 convol 3.01e- 2
## 7 35 convol 2.04e- 1
## 8 40 convol 8.84e- 1
## 9 45 convol 2.77e+ 0
## 10 50 convol 6.80e+ 0
## # ... with 110 more rowsggplot(mixed_table2, aes(x = age, y = volume, colour = type)) +
geom_line() +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)") +
ggtitle("Mixedwood yield curves")
Create a table of yields for pine/black spruce plots
age <- 1:40 * 5
plsb_table <- tibble(age, totvol = cr_volume(300.0,0.0275,7.068,age),
convol = cr_volume(300,0.0190,4.371,age),
decvol = totvol - convol)
plsb_table## # A tibble: 40 x 4
## age totvol convol decvol
## <dbl> <dbl> <dbl> <dbl>
## 1 5 0.000151 0.00830 -0.00815
## 2 10 0.0126 0.140 -0.128
## 3 15 0.140 0.676 -0.536
## 4 20 0.686 1.95 -1.27
## 5 25 2.15 4.27 -2.13
## 6 30 5.09 7.85 -2.75
## 7 35 10.0 12.8 -2.77
## 8 40 17.2 19.1 -1.90
## 9 45 26.6 26.6 -0.0214
## 10 50 38.2 35.4 2.84
## # ... with 30 more rowsplsb_table2 <-
plsb_table |>
pivot_longer(cols=c("totvol","convol","decvol"), names_to = "type",
values_to = "volume") |>
arrange(type,age)
plsb_table2## # A tibble: 120 x 3
## age type volume
## <dbl> <chr> <dbl>
## 1 5 convol 0.00830
## 2 10 convol 0.140
## 3 15 convol 0.676
## 4 20 convol 1.95
## 5 25 convol 4.27
## 6 30 convol 7.85
## 7 35 convol 12.8
## 8 40 convol 19.1
## 9 45 convol 26.6
## 10 50 convol 35.4
## # ... with 110 more rowsggplot(plsb_table2, aes(x = age, y = volume, colour = type)) +
geom_line() +
ylim(0,600) +
xlab("age (y)") +
ylab("volume (m³/ha)") +
ggtitle("Lodgepole pine/black spruce yield curves")## Warning: Removed 9 rows containing missing values (`geom_line()`).
Make some pretty tables
mixed_table |>
gt() |>
tab_header(title = "Yield tables for aspen/white spruce stands") |>
fmt_number(columns = c("totvol","convol","decvol"),decimals = 1)| Yield tables for aspen/white spruce stands | |||
| age | totvol | convol | decvol |
|---|---|---|---|
| 5 | 0.0 | 0.0 | 0.0 |
| 10 | 0.0 | 0.0 | 0.0 |
| 15 | 0.0 | 0.0 | 0.0 |
| 20 | 0.4 | 0.0 | 0.4 |
| 25 | 2.6 | 0.0 | 2.6 |
| 30 | 9.8 | 0.0 | 9.7 |
| 35 | 25.3 | 0.2 | 25.1 |
| 40 | 50.2 | 0.9 | 49.3 |
| 45 | 82.9 | 2.8 | 80.2 |
| 50 | 120.2 | 6.8 | 113.4 |
| 55 | 158.4 | 13.9 | 144.5 |
| 60 | 194.6 | 24.5 | 170.0 |
| 65 | 227.0 | 38.8 | 188.3 |
| 70 | 255.0 | 56.1 | 198.9 |
| 75 | 278.2 | 75.5 | 202.7 |
| 80 | 297.2 | 96.2 | 201.0 |
| 85 | 312.3 | 117.1 | 195.2 |
| 90 | 324.3 | 137.3 | 186.9 |
| 95 | 333.6 | 156.4 | 177.2 |
| 100 | 340.8 | 174.0 | 166.9 |
| 105 | 346.4 | 189.7 | 156.7 |
| 110 | 350.7 | 203.6 | 147.1 |
| 115 | 354.0 | 215.8 | 138.2 |
| 120 | 356.5 | 226.2 | 130.2 |
| 125 | 358.4 | 235.2 | 123.2 |
| 130 | 359.8 | 242.7 | 117.1 |
| 135 | 360.9 | 249.1 | 111.8 |
| 140 | 361.8 | 254.4 | 107.3 |
| 145 | 362.4 | 258.9 | 103.5 |
| 150 | 362.9 | 262.6 | 100.3 |
| 155 | 363.3 | 265.7 | 97.6 |
| 160 | 363.5 | 268.2 | 95.3 |
| 165 | 363.7 | 270.3 | 93.4 |
| 170 | 363.9 | 272.1 | 91.8 |
| 175 | 364.0 | 273.5 | 90.5 |
| 180 | 364.1 | 274.7 | 89.4 |
| 185 | 364.2 | 275.7 | 88.5 |
| 190 | 364.2 | 276.5 | 87.8 |
| 195 | 364.3 | 277.1 | 87.2 |
| 200 | 364.3 | 277.7 | 86.7 |
plsb_table |>
gt() |>
tab_header(title = "Yield tables for pine/black spruce stands") |>
fmt_number(columns = c("totvol","convol","decvol"),decimals = 1)| Yield tables for pine/black spruce stands | |||
| age | totvol | convol | decvol |
|---|---|---|---|
| 5 | 0.0 | 0.0 | 0.0 |
| 10 | 0.0 | 0.1 | −0.1 |
| 15 | 0.1 | 0.7 | −0.5 |
| 20 | 0.7 | 2.0 | −1.3 |
| 25 | 2.1 | 4.3 | −2.1 |
| 30 | 5.1 | 7.8 | −2.8 |
| 35 | 10.0 | 12.8 | −2.8 |
| 40 | 17.2 | 19.1 | −1.9 |
| 45 | 26.6 | 26.6 | 0.0 |
| 50 | 38.2 | 35.4 | 2.8 |
| 55 | 51.6 | 45.1 | 6.5 |
| 60 | 66.5 | 55.7 | 10.8 |
| 65 | 82.2 | 66.8 | 15.4 |
| 70 | 98.4 | 78.3 | 20.1 |
| 75 | 114.7 | 90.1 | 24.6 |
| 80 | 130.8 | 102.0 | 28.8 |
| 85 | 146.3 | 113.8 | 32.5 |
| 90 | 161.2 | 125.4 | 35.7 |
| 95 | 175.1 | 136.8 | 38.3 |
| 100 | 188.1 | 147.8 | 40.3 |
| 105 | 200.1 | 158.3 | 41.7 |
| 110 | 211.0 | 168.5 | 42.6 |
| 115 | 221.0 | 178.1 | 42.9 |
| 120 | 230.0 | 187.2 | 42.8 |
| 125 | 238.1 | 195.8 | 42.3 |
| 130 | 245.4 | 203.9 | 41.5 |
| 135 | 251.9 | 211.4 | 40.5 |
| 140 | 257.7 | 218.5 | 39.2 |
| 145 | 262.8 | 225.1 | 37.7 |
| 150 | 267.4 | 231.2 | 36.1 |
| 155 | 271.4 | 236.9 | 34.5 |
| 160 | 274.9 | 242.1 | 32.8 |
| 165 | 278.0 | 247.0 | 31.0 |
| 170 | 280.8 | 251.5 | 29.3 |
| 175 | 283.2 | 255.6 | 27.6 |
| 180 | 285.3 | 259.4 | 25.9 |
| 185 | 287.2 | 262.9 | 24.3 |
| 190 | 288.8 | 266.1 | 22.7 |
| 195 | 290.2 | 269.1 | 21.1 |
| 200 | 291.4 | 271.8 | 19.7 |