PCA covariate reduction for prediction grid

[kyle-messier]

We need to think through how PCA reduction of TRI variables will impact the prediction grid covariates. ie. We need to ensure that the same PCA variables can be calculated in the prediction locations.

Does that mean we retain the PCA loading from fit locations and apply to predictions or do separate/new PCA and retain the same # of PCs (the latter would presumably have different loadings - is that a big deal?)

[kyle-messier]

recipes should be able to handle this. We need to apply the SAME transformation to the prediction locations.

library(recipes)

rec <- recipe(y ~ ., data = training_data) %>%
  step_normalize(all_predictors()) %>%
  step_pca(all_predictors(), num_comp = 5)

rec_prep <- prep(rec, training = training_data)

X_train_pca <- bake(rec_prep, new_data = training_data)
X_pred_pca  <- bake(rec_prep, new_data = prediction_data)

[MAKassien]

Another option could be collinearity reduction based on correlation between different TRI columns. This helps cut down on the number of TRI covariates but keeps them more traceable for reproducibility in the predictions.

I did a preliminary check on this with 2021 yearly TRI data from my HAPSGP project, and when cutting collinearity using this function I built for my project, the amount of TRI columns were reduced from 2125 to 355. This is a significant cut although it is still way more than what we can achieve using PCA (currently we only use 10 TRI columns, from 5 principal components at two different radii)

[mitchellmanware]

@kyle-messier @MAKassien
TRI PCA reduction uses beethoven::reduce_pca function, which "bakes" a "prepped" recipe to new data. If we save this prepared recipe in the dt_feat_calc_tri target where it is run, we should be able to "bake" the prediction grid TRI covariates in the same way.

beethoven/R/calc_postprocessing.R

Line 822 in e9327ba

data_prep <- recipes::prep(data_pca, data = data)

Updating the dt_feat_calc_tri target will invalidate the imputed covariate data, so I will work on this change once we know the base/meta learner stage updates are working.

[mitchellmanware]

@kyle-messier @sigmafelix
Review of the reduce_pca and post_calc_pca for extracting the recipe has shown that the groups parameter of post_calc_pca sometimes catches TRI columns outside of the radius "group". Current setup is designed to perform PCA on predictors present within the 1000, 10000, 50000 meter range, but using a simple grep("01000", names(dt_tri)) sometimes includes 10000 or 50000 if "01000" is part of the dynamic branch's random name. Example below should be only 1000 radius predictors.

[377] "STACK_AIR_0000606202_01000"    "STACK_AIR_0000075025_01000"   
[379] "STACK_AIR_0000081072_01000"    "FUGITIVE_AIR_0000100027_10000"
[381] "STACK_AIR_0000100027_10000"    "STACK_AIR_0000100016_10000"   
[383] "FUGITIVE_AIR_0000100027_50000" "FUGITIVE_AIR_0000100016_50000"
[385] "STACK_AIR_0000100027_50000"    "STACK_AIR_0000100016_50000"   
[387] "FUGITIVE_AIR_N535_01000"       "STACK_AIR_N535_01000"         
> 

I will correct this with a more robust regular expression, but will require re-dispatching dt_feat_calc_triand subsequent targets. We needed to do it anyway to extract the recipe, so it is not terrible. Also, since they are latent variables and only 15/300+ included in the base/meta learners, the results of the models should not change much - we can still interpret the current results and adjust as necessary.

[mitchellmanware]

Retention of recipe and post-bake returns same values, so application to the prediction grid tri values should be good with some updates for formatting column names for the 01000 10000 50000 groups.

> dt_pca2 <- beethoven::post_calc_pca(
+   data = dt_wide,
+   locs_id = "site_id",
+   time_id = "time",
+   yvar = "Arithmetic.Mean",
+   coords = c("lon", "lat"),
+   num_comp = 5,
+   threshold = NA,
+   pattern = "FUGITIVE|STACK",
+   groups = c("01000", "10000", "50000"),
+   prefix = "TRI"
+ )
> dt_pca_post <- recipes::bake(dt_pca2$recipe[[1]], new_data = dt_wide)
> all.equal(
+   unname(data.table::data.table(dt_pca_post)),
+     dt_pca2$pca[, grep("01000", names(dt_pca2$pca), value = TRUE), with = FALSE]
+   )
+ )
[1] TRUE

[mitchellmanware]

be40820
Application to dt_tri_feat_calc incoming.