# Accounting for measurement error during variable selection with projpred (possibly with rstanarm or brms?)

**URL:** <https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789>\
**Category:** Modeling\
**Created:** [November 7, 2019, 3:32am UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789 "2019-11-07T03:32:43Z")\
**Posts on this page:** 6\
**Page:** 2

<div class="post-metadata">

**Author:** ![fweber144](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/fweber144/32/15592_2.png) [@fweber144](https://discourse.mc-stan.org/u/fweber144)\
**Post date:** [June 24, 2025, 10:54am UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/21 "2025-06-24T10:54:28Z")

</div>

> [@rtnliqry](#):
>
> This now works when `parallel = FALSE` but throws the following error when parallelised:
> 
> > Error in { : task 1 failed - “could not find function “nlist””

You might be able to resolve this by using the most recent GitHub version and setting `options(projpred.export_to_workers = c("loo::nlist"))` or (if that does not work) assigning `nlist_glob <- loo::nlist` and then setting `options(projpred.export_to_workers = c("nlist_glob"))`.

> [@rtnliqry](#):
>
> When trying `latent = TRUE`, it unfortunately throws the same error

I’ll check that asap.

---

<div class="post-metadata">

**Author:** ![rtnliqry](https://avatars.discourse-cdn.com/v4/letter/r/df705f/32.png) [@rtnliqry](https://discourse.mc-stan.org/u/rtnliqry)\
**Post date:** [June 24, 2025, 2:23pm UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/22 "2025-06-24T14:23:49Z")

</div>

> [@fweber144](#):
>
> You might be able to resolve this by using the most recent GitHub version and setting `options(projpred.export_to_workers = c("loo::nlist"))`

This works, thanks!

> [@fweber144](#):
>
> Btw, the fix also shows that the ((I)S)PCA approach introduces a dependency between training and test data: When making predictions for the left-out observations from fold `k`, we use the ((I)S)PCA `predict()` method which relies on estimates based on the training data for fold `k`

Is this not still making predictions on data neither the (I)SPCA nor the `brms` model has seen before? **Edit** : Related to this, is there a difference in the way `varsel` handles this compared to `cv_varsel`, given that it works even if you don’t supply a custom `ref_predfun`?

> [@fweber144](#):
>
> And of course, _supervised_ PCA (the iterative variant as well) also comes with a dependency between response and predictors, but I think that’s clear already from the name (and this has been pointed out by Piironen and Vehtari (2018): …

Yes, I’ve spent a long time thinking about this, but it does seem equivalent to a standard machine learning workflow, where one would (assuming no hyperparameters are being tuned): (i) select features within each fold of a K-fold CV loop; (ii) build a model for each fold using the selected feature subset (on the same data); and then (iii) predict on the left-out data. I assume as long as one predicts on data neither the feature selection process nor model have ever seen, the K-fold CV performance estimate should be unbiased?

> [@fweber144](#):
>
> And there might also be other relevant papers apart from Piironen and Vehtari (2018).

The out-of-sample performance for feature selection is explored in [Pavone _et al_ (2023)](https://link.springer.com/article/10.1007/s00180-022-01231-6) for the `bodyfat` data (code [here](https://github.com/fpavone/ref-approach-paper/blob/master/code/bodyfat_comparison_outofsample.R)). They seem to nest an exact LOO-CV within a bootstrap to estimate the predicted outcome value (`yfit`) for each observation using a reference model constructed using all the other data. I’m not sure how this could be done within `projpred` – I suppose one could perform an inner K-fold CV loop within each iteration of `cvfun`?

I’d also welcome @avehtari’s thoughts on this!

---

<div class="post-metadata">

**Author:** ![fweber144](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/fweber144/32/15592_2.png) [@fweber144](https://discourse.mc-stan.org/u/fweber144)\
**Post date:** [June 25, 2025, 10:19am UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/23 "2025-06-25T10:19:40Z")

</div>

> [@fweber144](#):
>
> > [@rtnliqry](#):
> >
> > When trying `latent = TRUE`, it unfortunately throws the same error
> 
> I’ll check that asap.

This was indeed a deficiency in the internal default `extract_model_data` function which is defined in `init_refmodel()`. I fixed this in the current GitHub version (branch `master`, see also PR [#523](https://github.com/stan-dev/projpred/pull/523)). Thanks for the report and the good reprex!

---

<div class="post-metadata">

**Author:** ![fweber144](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/fweber144/32/15592_2.png) [@fweber144](https://discourse.mc-stan.org/u/fweber144)\
**Post date:** [July 3, 2025, 7:42am UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/24 "2025-07-03T07:42:38Z")

</div>

> [@rtnliqry](#):
>
> Is this not still making predictions on data neither the (I)SPCA nor the `brms` model has seen before?

It depends a bit on what you consider as data: (i) only response values, (ii) response and original predictor values, or (iii) response and ((I)S)PC-constructed predictor values.

Above, I had (iii) in mind when I mentioned the dependency between training and test data. For (iii), the brms model has indeed seen parts of this data before because the ((I)S)PC-constructed predictor values for the left-out observations include information from the training data. However, now that you are asking this, I’m realizing that this dependency might not be an issue because the theory behind projpred _conditions_ on the predictor data (just like many other types of analyses do). So it was probably not correct from my side to think of (iii) as “data” in this case (although I didn’t do this on purpose; I was simply just thinking of “data” as the dataset and didn’t think it through).

For (i) and (ii), you would be right that this data has not been seen by the trained model when making predictions for the left-out observations.

> [@rtnliqry](#):
>
> **Edit** : Related to this, is there a difference in the way `varsel` handles this compared to `cv_varsel`, given that it works even if you don’t supply a custom `ref_predfun`?

Yes, there is a difference which explains that: Since `varsel()` only performs in-sample predictions, we never get `!is.null(newdata)` in `ref_predfun`, so the adjustment I suggested above for K-fold CV is irrelevant. Hence the internal default `ref_predfun` (which is very similar to your `custom_ref_predfun`) does not throw an error in the `varsel()` case.

> [@rtnliqry](#):
>
> I assume as long as one predicts on data neither the feature selection process nor model have ever seen, the K-fold CV performance estimate should be unbiased?

Yes, for a given submodel size, the K-fold CV performance estimate should be unbiased in that case.

> [@rtnliqry](#):
>
> They seem to nest an exact LOO-CV within a bootstrap to estimate the predicted outcome value (`yfit`) for each observation using a reference model constructed using all the other data. I’m not sure how this could be done within `projpred` – I suppose one could perform an inner K-fold CV loop within each iteration of `cvfun`?

I don’t think this is necessary for typical projpred use-cases. In principle, an inner CV nested in an outer CV could help against the problem that the selection of the final submodel size is not cross-validated (for this, the selection of the final submodel size would have to be performed in each iteration of the outer CV loop—which is not straightforward if the size selection is a manual one), but I’m not sure how a CV nested within a bootstrap should help. Such a CV nested within a bootstrap could probably be used to assess the sampling distribution of a predictive performance (e.g., ELPD) estimator (at each submodel size separately), but I don’t think this is necessary in typical projpred use-cases. It has been a while since I read Pavone et al. (2023), so I currently cannot say more about this, sorry.

---

<div class="post-metadata">

**Author:** ![rtnliqry](https://avatars.discourse-cdn.com/v4/letter/r/df705f/32.png) [@rtnliqry](https://discourse.mc-stan.org/u/rtnliqry)\
**Post date:** [July 29, 2025, 12:38pm UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/25 "2025-07-29T12:38:06Z")

</div>

Thanks for your response above, @fweber144, much appreciated.

Somewhat tangentially related to this: say a `cvfits` list is a huge object because the individual model fits are large, does it make sense to run multiple individual `cv_varsel` runs, each with the same reference model fit but with a different subset of the `cvfits` list, and then combine the results afterwards?

---

<div class="post-metadata">

**Author:** ![fweber144](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/fweber144/32/15592_2.png) [@fweber144](https://discourse.mc-stan.org/u/fweber144)\
**Post date:** [July 31, 2025, 8:45pm UTC](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789/26 "2025-07-31T20:45:01Z")

</div>

I’m not sure if this is possible in projpred, but feel free to try it out.

[Previous page](https://discourse.mc-stan.org/t/accounting-for-measurement-error-during-variable-selection-with-projpred-possibly-with-rstanarm-or-brms/11789.md?page=1)