Generating all possible interaction terms for model building

Question

This was the question which was solved by fellow stackoverflow users.

To build the regression model for each PCs and the independent predictors this was used

iv <- c("Sex", "fAge", "Index", "Lane", "Gen")  
dv <- paste0('PC', 1:8)
rhs <- unlist(sapply(1:length(iv), function(m) apply(combn(iv, m = m), 2, paste, collapse = ' + ')))
frms <- with(expand.grid(dv, rhs), paste(Var1, Var2, sep = ' ~ '))

models <- lapply(frms, function(x) anova(lm(x, data = mrna.pcs)))
names(models) <- frms

Few of the models out which i ran in to test are as such

"PC3 ~ Sex + fAge"                      "PC4 ~ Sex + fAge"                      "PC5 ~ Sex + fAge"  

I was read through how to use interaction terms from one of the post so I realised I'm missing something like this

y ~ r*s
r + s + rs

In my case i would like to see something this sort as one of the example.

PC6 ~ Sex + fAge + Sex:fAge 

So how do I generate the all possible interaction term from predictors

Now what I have is no interaction terms in my models.

Any suggestion or help would be really appreciated. Answer Updated

@Jay.sf suggested example

lm(mpg ~ (cyl+disp+hp)^2, data = mtcars)

This works exactly what I would expect to be

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  5.601e+01  8.556e+00   6.546  7.4e-07 ***
cyl         -4.427e+00  2.474e+00  -1.789   0.0857 .  
disp        -1.184e-01  4.780e-02  -2.476   0.0204 *  
hp          -1.142e-01  7.394e-02  -1.544   0.1352    
cyl:disp     1.439e-02  7.153e-03   2.012   0.0551 .  
cyl:hp       1.556e-02  1.364e-02   1.141   0.2647    
disp:hp     -8.567e-05  2.201e-04  -0.389   0.7005 

but when i use this as he suggested

rhs <- unlist(sapply(1:length(iv), function(m) apply(combn(iv, m = m), 2, paste, collapse = ' * ')))

The model output are as such

$`PC7 ~ Sex * fAge`

Call:
lm(formula = x, data = mrna.pcs)

Coefficients:
(Intercept)         SexM       fAge12       fAge18  SexM:fAge12  SexM:fAge18  
     2.2140      -0.8965      -2.8590      -5.1501       1.9537       3.2378  


$`PC8 ~ Sex * fAge`

Call:
lm(formula = x, data = mrna.pcs)

Coefficients:
(Intercept)         SexM       fAge12       fAge18  SexM:fAge12  SexM:fAge18  
     1.0590       0.9449      -2.0983      -0.5103      -1.2526      -2.5817  

Not as the one which coming from the example i ran. How do i incorporate both of them + and : in the same code that would fix the problem i think

I think i was running the wrong one here is my outout

Response: PC7
           Df  Sum Sq Mean Sq F value Pr(>F)
Sex         1    21.9  21.898  0.2065 0.6501
fAge        2   388.7 194.371  1.8329 0.1629
Sex:fAge    2    83.4  41.693  0.3932 0.6755
Residuals 182 19300.3 106.046               

$`PC8 ~ Sex * fAge`
Analysis of Variance Table

Response: PC8
           Df  Sum Sq Mean Sq F value Pr(>F)
Sex         1     7.4   7.393  0.0787 0.7794
fAge        2   236.5 118.261  1.2581 0.2866
Sex:fAge    2    52.3  26.138  0.2781 0.7576
Residuals 182 17107.4  93.996