How do I make the p_adjust function from hdm package work?

24 Views Asked by At

I’m trying to fit lasso regressions in R, using the rlassoEffects function from the hdm package. I'm new to lasso procedure but so far I have managed to make it work. However, I now would like to use the function p_adjust to adjust my p-values (with Romano-Wolf method), and I'm struggling with that.

I've managed to make it work once with my first lasso model, but not with a second one. I've received the following error message: "Error in colMeans(v^2) : 'x' must be an array of at least two dimensions". I guess the structure of my second model must be a bit different from the expected structure for input to the p_adjust function, but I don't see what's wrong.

This is what my first model looks like:

lasso_1 = hdm::rlassoEffects(x = covariates_matrix, y = outcome_vector, index = c(1:2), method = "double selection")

summary(lasso_1)

[1] "Estimates and significance testing of the effect of target variables"
             Estimate. Std. Error t value Pr(>|t|)
treatment1     0.11926    0.09154   1.303    0.193
treatment2     0.10241    0.08389   1.221    0.222

Then I ran this code p_adjust(lasso_1, method = "RW", B = 1000) and it worked fine.

And this is the second model:

lasso_2 = hdm::rlassoEffect(x = covariates_matrix, y = outcome_vector, d = treatment_vector, method = "double selection")

summary(lasso_2)

[1] "Estimates and significance testing of the effect of target variables"
   Estimate. Std. Error t value Pr(>|t|)  
d1   -0.3906     0.1965  -1.988   0.0468 *

And this time I got the error message when using p_adjust. I've tried to manually add an intercept term to my lasso_2 object, just to see if p_adjust works better once my model has a second row:


coefficients <- lasso_2$coefficients
se <- lasso_2$se
pval <- lasso_2$pval
t <- lasso_2$t

intercept <- c(intercept = 0)

# Combine intercept and coefficients
new_coefficients <- c(intercept, coefficients)
new_se <- c(intercept, se)
new_pval <- c(intercept, pval)
new_t <- c(intercept, t)


# Create a modified rlassoEffect object
lasso_2_with_intercept <- lasso_2

lasso_2_with_intercept$coefficients <- new_coefficients
lasso_2_with_intercept$se <- new_se
lasso_2_with_intercept$pval <- new_pval
lasso_2_with_intercept$t <- new_t

But when I run p_adjust(lasso_2_with_intercept, method = "RW", B = 1000), I still get the same error message... So it seems that p_adjust expects a specific data structure that's not provided in my lasso_2 object, even with a few quick modifications.

Do you think the error comes from differences in the specification of the lasso models? Or is it just because I have only one treatment variable in my second model? But surely p_adjust also works when you just have one treatment variable? Or maybe some of you can recommend another package that works better to adjust the p-values of lasso regressions?

I'm not including a reproducible example because my code is very long and based on a very large dataset, but I'm hoping that someone can detect some obvious mistake based on what I've mentioned, because I'm quite new to these lasso procedures and not a pro in R.

0

There are 0 best solutions below