Checking out SLA and Vcmax

Objective

1. Compare default ED with UMBS SLA and Vcmax
2. Update ED SLA to UMBS values

Markdown set up

1. Ed vs UMBS

Compare the default ED SLA and Vcmax values with observations from UMBS.

A. Import Data

Import the data sets.

# The results from  a ED run 
ed_data <- readRDS(file.path(OUTPUT_DIR, "full-length.rds"))
# Read UMBS data from the fortedata  pacakge. 
umbs_sla <- as.data.table(read.csv(system.file('extdata/fd_sla.csv', package = 'fortedata'), 
                                   stringsAsFactors = FALSE))

Create mapping data frame that will be used to map ed pfts to the umbs species.

# The ED pfts 
PFTs    <- c(6, 8, 9, 10, 11, 7) # The PFTS to keep, some of the ED outputs report all of the PFTs.
pft_names <- data.table(pft = PFTs, 
                        pft_name = c('Northern North American pines', 
                                     'Late-successional North American conifers', 
                                     'Temperate broadleaf, early successional',  
                                     'Temperate broadleaf, mid-successional', 
                                     'Temperate broadleaf, late successional', 
                                     'Southern North American pines'), stringsAsFactors = FALSE)
pft_names <- pft_names[! pft %in% c(7, 8)]
# UMBS to ED pfts
mid_success       <- data.table(pft = 10, UMBS = c('ACRU','ACPE',  'FAGR', 'ACSA'), stringsAsFactors = FALSE)
temp_broadl_late  <- data.table(pft = 11, UMBS = c('QURU', 'TCSA'), stringsAsFactors = FALSE)
nor_pines         <- data.table(pft = 6, UMBS = c('PIST', 'PIRE'), stringsAsFactors = FALSE)
temp_broadl_early <- data.table(pft = 9, UMBS = c("AMEL", "POGR", "BEPA", "BEAL", "POTR"), stringsAsFactors = FALSE)
umbs_ed <- rbind(mid_success, temp_broadl_late, nor_pines, temp_broadl_early)

Let’s take a look at the maping betwee the UMBS species and the ED pfts.

umbs_ed[pft_names, on = 'pft', nomatch = 0][ , list(pft_name, pft, UMBS)] %>% 
  knitr::kable(format = 'markdown')

pft_name	pft	UMBS
Northern North American pines	6	PIST
Northern North American pines	6	PIRE
Temperate broadleaf, early successional	9	AMEL
Temperate broadleaf, early successional	9	POGR
Temperate broadleaf, early successional	9	BEPA
Temperate broadleaf, early successional	9	BEAL
Temperate broadleaf, early successional	9	POTR
Temperate broadleaf, mid-successional	10	ACRU
Temperate broadleaf, mid-successional	10	ACPE
Temperate broadleaf, mid-successional	10	FAGR
Temperate broadleaf, mid-successional	10	ACSA
Temperate broadleaf, late successional	11	QURU
Temperate broadleaf, late successional	11	TCSA

B. Format UMBS SLA

# Add the ED PFT information to the UMBS SLA measurements. 
umbs_sla <- umbs_sla[umbs_ed, on = c('species' = 'UMBS'), nomatch  = 0]

As Alexey noted, convert the SLA units.

umbs_sla$sla <- umbs_sla$sla * (1/0.48)

# Now determine the average value for each PFT, this is more analagous to what ED  does. 
umbs_sla <- umbs_sla[ , list(sla = mean(sla), min = min(sla), max = max(sla)), by = list(pft)]
umbs_sla <- pft_names[umbs_sla, on = "pft", nomatch = 0][, list(pft_name, pft, sla, min, max)] 

Let’s take a look at what UMBS SLA values look like

knitr::kable(umbs_sla)

pft_name	pft	sla	min	max
Temperate broadleaf, mid-successional	10	47.91667	39.58333	72.91667
Temperate broadleaf, late successional	11	20.87500	12.16667	29.58333
Northern North American pines	6	18.63542	11.22917	26.04167
Temperate broadleaf, early successional	9	37.46667	33.10417	43.37500

C. Format ED SLA

Import ED results and select the SLA values from our baseline run. Note that we will only be able to extract SLA values from the output, the Vcmax will have to come in put parameters.

ed_data$df_cohort[[1]] %>%
  dplyr::select(SLA, datetime, DBH, PFT) %>%  
  dplyr::left_join(pft_names, by = c('PFT' = "pft")) %>%  
  dplyr::select(pft_name, pft = PFT, SLA)  %>% 
  dplyr::distinct() %>%  
    as.data.table() -> 
    ed_sla
  knitr::kable(ed_sla)

pft_name	pft	SLA
Northern North American pines	6	6.0
Temperate broadleaf, mid-successional	10	24.2
Temperate broadleaf, late successional	11	60.0
Temperate broadleaf, early successional	9	30.0

D. Compare ED and UMBS SLA

ggplot() + 
    geom_point(data = umbs_sla, aes(pft_name, min, color = 'umbs', shape = 'min sla'),  size = 4) + 
    geom_point(data = umbs_sla, aes(pft_name, max, color = 'umbs', shape = 'max sla'),  size = 4) + 
    geom_point(data = umbs_sla, aes(pft_name, sla, color = 'umbs', shape = 'mean sla'),  size = 4) + 
  geom_point(data = ed_sla, aes(pft_name, SLA, color = 'ed', shape = 'mean sla'),  size = 4) + 
  THEME + 
  theme( axis.text.x = element_text(angle = 25, hjust = 1), 
         legend.position = 'bottom', 
         legend.text = element_text(size = 16)) + 
  labs(title = 'Comparison of UMBS and ED SLA', 
       y = 'SLA', 
       x = 'ED PFTs')

Recall from baseline ed run that the stand is dominated by the temperate broadlead early succesional.

E. ED Vcmax

Taking a look at ED code we can pull out the \({V}_{m0}\) values which is \({V}_{cmax}\) at 15 C instead of 25 C see Alexey’s notes.

ed_Vm0 <- cbind(pft_names, Vm0 = c(11.350000, 20.387075, 17.454687, 6.981875))

Now need to convert from the 15C to 25C, we can use the arrhenius scaling function from the fortebaseline.

ed_Vcmax_values <- arrhenius.scaling(observed.value = ed_Vm0$Vm0, 
                                     old.temp = 15,
                                     new.temp = 25)
ed_Vcmax <- cbind(pft_names, Vcmax = ed_Vcmax_values)

F. UMBS Vcmax

Lisa sent me some data for UMBS to provide a ball park estimate of the Vcmax values.

UMBS data table

UMBS  <- c('ACRU', 'BEPA', 'QURU', 'POGR')
Vcmax <- c(47.3, 58.5, 62.6, 71.4)
umbs_vcmax <- data.table(UMBS = UMBS, Vcmax = Vcmax)
umbs_vcmax <- umbs_vcmax[umbs_ed, on = 'UMBS', nomatch = 0]
umbs_vcmax <- umbs_vcmax[pft_names, on = 'pft', nomatch = 0]

G. Vcmax Comparison

Compare the UMBS and the ED Vcmax values.

ggplot() + 
  geom_point(data = ed_Vcmax, aes(pft_name, Vcmax, color = 'ED'), size = 3) + 
    geom_point(data = umbs_vcmax, aes(pft_name, Vcmax, color = 'UMBS'), size = 3) + 
    THEME + 
  theme( axis.text.x = element_text(angle = 25, hjust = 1), 
         legend.position = 'bottom', 
         legend.text = element_text(size = 16)) + 
  labs(title = 'Comparison of UMBS and ED Vcmax', 
       y = 'Vcmax', 
       x = 'ED PFTs')

Vcmax seems to be fairly high.

H. Concluding Thoughts & Next Steps

• This analysis only makes sense if the ED PFTs and UMBS species are categorized properly
• The ED temprate borad leaf late sucessional SLA is really far off! Could this be throwing the ED forest off?

1. Ed run with “realistic” SLA
2. Ed run with “realistic” Vcmax
3. Ed run with “realistic” SLA and Vcmax

2. Adjust ED SLA

What happens when the ED SLA values are updated to better reflect the UMBS SLA values?

A. Import and format the results from the SLA adjusted run.

Start by importing and formating results from the FoRTE baseline run.

# Read the input file 
baseline <- readRDS(file.path(OUTPUT_DIR, 'test', 'baseline.rds'))
# Extract the SLA values, these will be compared with the old sla ED values and the UMBS observations. 
baseline$df_cohort[[1]] %>%
  dplyr::select(SLA, datetime, DBH, PFT) %>%  
  dplyr::left_join(pft_names, by = c('PFT' = "pft")) %>%  
  dplyr::select(pft_name, pft = PFT, SLA)  %>% 
  dplyr::distinct() %>%  
  as.data.table() ->
 ED_default_sla
ED_default_sla$scn <- 'default SLA'
ED_default_LAI <- as.data.table(baseline$df_pft[[1]])[  ,list(datetime, MMEAN_LAI_PY, pft)]
ED_default_LAI$datetime <- ymd(ED_default_LAI$datetime)
ED_default_LAI$year     <- year(ED_default_LAI$datetime)
ED_default_LAI$month    <- month(ED_default_LAI$datetime)
ED_default_LAI$variable <- 'MMEAN_LAI_CO'
ED_default_LAI <- ED_default_LAI[as.data.table(ed2_variable_info())[, list(variable, description, unit)], 
           on = 'variable',  nomatch=0]
ED_default_LAI     <- ED_default_LAI[pft_names, on = c('pft')]
ED_default_LAI$scn <- 'default SLA'
ED_default_NPP <- as.data.table(baseline$df_scalar[[1]])[  ,list(datetime, MMEAN_NEP_PY)]
ED_default_NPP$datetime <- ymd(ED_default_NPP$datetime)
ED_default_NPP$year <- year(ED_default_NPP$datetime)
ED_default_NPP$month <- month(ED_default_NPP$datetime)
ED_default_NPP$variable <- 'MMEAN_NEP_PY'
ED_default_NPP <- ED_default_NPP[as.data.table(ed2_variable_info())[, list(variable, description, unit)], 
                                 on = 'variable',  nomatch=0]
ED_default_NPP <- ED_default_NPP[, list(value = mean(MMEAN_NEP_PY, na.rm = TRUE)), 
                                 by = list(year, variable, unit, description)]
ED_default_NPP$scn <- 'ED default'

B. Import and format the results from the SLA adjusted run.

# Read in the input 
umbs_sla_adjusted <- readRDS(file.path(OUTPUT_DIR, 'test', 'UMBS_SLA-1.rds'))
# Extract the SLA values, these will be compared with the old sla ED values and the UMBS observations. 
umbs_sla_adjusted$df_cohort[[1]] %>%
  dplyr::select(SLA, datetime, DBH, PFT) %>%  
  dplyr::left_join(pft_names, by = c('PFT' = "pft")) %>%  
  dplyr::select(pft_name, pft = PFT, SLA)  %>% 
  dplyr::distinct() %>%  
  as.data.table() ->
  ED_adjusted_sla
ED_adjusted_sla$scn <- 'adjusted SLA'
umbs_sla_adjusted_LAI <- as.data.table(umbs_sla_adjusted$df_pft[[1]])[  ,list(datetime, MMEAN_LAI_PY, pft)]
umbs_sla_adjusted_LAI$datetime <- ymd(umbs_sla_adjusted_LAI$datetime)
umbs_sla_adjusted_LAI$year     <- year(umbs_sla_adjusted_LAI$datetime)
umbs_sla_adjusted_LAI$month    <- month(umbs_sla_adjusted_LAI$datetime)
umbs_sla_adjusted_LAI$variable <- 'MMEAN_LAI_CO'
umbs_sla_adjusted_LAI <- umbs_sla_adjusted_LAI[as.data.table(ed2_variable_info())[, list(variable, description, unit)], 
                                               on = 'variable',  nomatch=0]
umbs_sla_adjusted_LAI <- umbs_sla_adjusted_LAI[pft_names, on = c('pft')]
umbs_sla_adjusted_LAI$scn <- 'adjusted SLA'
umbs_ed_NPP <- as.data.table(umbs_sla_adjusted$df_scalar[[1]])[  ,list(datetime, MMEAN_NEP_PY)]
umbs_ed_NPP$datetime <- ymd(umbs_ed_NPP$datetime)
umbs_ed_NPP$year <- year(umbs_ed_NPP$datetime)
umbs_ed_NPP$month <- month(umbs_ed_NPP$datetime)
umbs_ed_NPP$variable <- 'MMEAN_NEP_PY'
umbs_ed_NPP <- umbs_ed_NPP[as.data.table(ed2_variable_info())[, list(variable, description, unit)], 
                           on = 'variable',  nomatch=0]
umbs_ed_NPP <- umbs_ed_NPP[, list(value = mean(MMEAN_NEP_PY, na.rm = TRUE)), 
                           by = list(year, variable, unit, description)]
umbs_ed_NPP$scn <- 'UMBS ed'
ED_umbs_AGB <- as.data.table(umbs_sla_adjusted$df_cohort[[1]])[  ,list(datetime, AGB_CO, PFT, DBH, NPLANT)]
# Add a unique identifier for each of the cohorts (reminder that cohorts = unique species at DBH)
split(ED_umbs_AGB, ED_umbs_AGB$datetime) %>%  
  lapply(function(x){
    
    x$CO <- LETTERS[1:nrow(x)]
    return(x)
    
  }) %>% 
  rbindlist() -> 
  ED_umbs_AGB
ED_umbs_AGB <- ED_umbs_AGB[pft_names, on = c('PFT' = 'pft')]
ED_umbs_AGB$year <- year(ED_umbs_AGB$datetime)
ED_umbs_AGB$month <- month(ED_umbs_AGB$datetime)
ED_umbs_AGB$variable <- "AGB_CO"
ED_umbs_AGB <- ED_umbs_AGB[as.data.table(ed2_variable_info())[, list(variable, description, unit)], on = 'variable',  nomatch=0]
ED_umbs_AGB$value <- ED_umbs_AGB$NPLANT * ED_umbs_AGB$AGB_CO
ED_umbs_AGB$unit <- gsub(pattern = '/plant', replacement = '/m2', ED_umbs_AGB$unit)
ED_umbs_AGB$scn <- 'umbs sla ed'

C. Compare Ed

The SLA values

ggplot() + 
  geom_point(data = umbs_sla, aes(pft_name, min, color = 'umbs', shape = 'min sla'),  size = 4) + 
    geom_point(data = umbs_sla, aes(pft_name, max, color = 'umbs', shape = 'max sla'),  size = 4) + 
    geom_point(data = umbs_sla, aes(pft_name, sla, color = 'umbs', shape = 'mean sla'),  size = 4) +
  geom_jitter(data = ED_adjusted_sla, aes(pft_name, SLA, color = 'umbs ED'), height = NULL,
              size = 4) + 
  geom_point(data = ED_default_sla, aes(pft_name, SLA, color = 'default ED'), size = 4) + 
  THEME + 
  theme( axis.text.x = element_text(angle = 25, hjust = 1), 
         legend.position = 'bottom', 
         legend.text = element_text(size = 16)) + 
  labs(title = 'Comparison of UMBS and ED SLA', 
       y = 'SLA', 
       x = 'ED PFTs')

Alright so the umbs ED SLA values are what we would expect! Which is awesome.

Biomass (above ground)

So it looks like by increasing the the SLA increased the biomass, expect for with the late sucessional temperate broadleaf which makes sense because we decreased the default SLA. The pines increased a lot…

NPP values

Hmmm it is intersting that the long term NPP is similar in the two scenarios, I guess that makes sense.

LAI values

Interesting the default LAI is larger in the long term, for the temprate broad leaf. While the LAI for the pines is substaintally larger.

D. Concluding Thoughts

So it surprises me that the adjusted SLA does not have a larger LAI long term for broadlead PFTs.