Search code examples
rdata.tablelinear-regression

Linear regression on dynamic groups in R


I have a data.table data_dt on which I want to run linear regression so that user can choose the number of columns in groups G1 and G2 using variable n_col. The following code works perfectly but it is slow due to extra time spent on creating matrices. To improve the performance of the code below, is there a way to remove Steps 1, 2, and 3 altogether by tweaking the formula of lm function and still get the same results?

library(timeSeries)
library(data.table)
data_dt = as.data.table(LPP2005REC[, -1])
n_col = 3 # Choose a number from 1 to 3
######### Step 1 ######### Create independent variable
xx <- as.matrix(data_dt[, "SPI"]) 
######### Step 2 ######### Create Group 1 of dependent variables
G1 <- as.matrix(data_dt[, .SD, .SDcols=c(1:n_col + 2)]) 
######### Step 3 ######### Create Group 2 of dependent variables
G2 <- as.matrix(data_dt[, .SD, .SDcols=c(1:n_col + 2 + n_col)]) 
lm(xx ~ G1 + G2)

Results -

summary(lm(xx ~ G1 + G2))
Call:
lm(formula = xx ~ G1 + G2)

Residuals:
       Min         1Q     Median         3Q        Max 
-3.763e-07 -4.130e-09  3.000e-09  9.840e-09  4.401e-07 

Coefficients:
              Estimate Std. Error    t value Pr(>|t|)    
(Intercept) -4.931e-09  3.038e-09 -1.623e+00   0.1054    
G1LMI       -5.000e-01  4.083e-06 -1.225e+05   <2e-16 ***
G1MPI       -2.000e+00  4.014e-06 -4.982e+05   <2e-16 ***
G1ALT       -1.500e+00  5.556e-06 -2.700e+05   <2e-16 ***
G2LPP25      3.071e-04  1.407e-04  2.184e+00   0.0296 *  
G2LPP40     -5.001e+00  2.360e-04 -2.119e+04   <2e-16 ***
G2LPP60      1.000e+01  8.704e-05  1.149e+05   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 5.762e-08 on 370 degrees of freedom
Multiple R-squared:      1, Adjusted R-squared:      1 
F-statistic: 1.104e+12 on 6 and 370 DF,  p-value: < 2.2e-16

Solution

  • This may be easier by just creating the formula with reformulate

    out <- lm(reformulate(names(data_dt)[c(1:n_col + 2, 1:n_col + 2 + n_col)], 
         response = 'SPI'), data = data_dt)
    

    -checking

    > summary(out)
    
    Call:
    lm(formula = reformulate(names(data_dt)[c(1:n_col + 2, 1:n_col + 
        2 + n_col)], response = "SPI"), data = data_dt)
    
    Residuals:
           Min         1Q     Median         3Q        Max 
    -3.763e-07 -4.130e-09  3.000e-09  9.840e-09  4.401e-07 
    
    Coefficients:
                  Estimate Std. Error    t value Pr(>|t|)    
    (Intercept) -4.931e-09  3.038e-09 -1.623e+00   0.1054    
    LMI         -5.000e-01  4.083e-06 -1.225e+05   <2e-16 ***
    MPI         -2.000e+00  4.014e-06 -4.982e+05   <2e-16 ***
    ALT         -1.500e+00  5.556e-06 -2.700e+05   <2e-16 ***
    LPP25        3.071e-04  1.407e-04  2.184e+00   0.0296 *  
    LPP40       -5.001e+00  2.360e-04 -2.119e+04   <2e-16 ***
    LPP60        1.000e+01  8.704e-05  1.149e+05   <2e-16 ***
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    
    Residual standard error: 5.762e-08 on 370 degrees of freedom
    Multiple R-squared:      1, Adjusted R-squared:      1 
    F-statistic: 1.104e+12 on 6 and 370 DF,  p-value: < 2.2e-16