Can I use DRIMseq to model interaction effect of variables?
One way to see the interaction effect would be including it in the design matrix as explained in the limma tutorial. Something on the line of:

design_full &lt;- model.matrix(~ treatment*strain, data = samples(d))


Or you could build it yourself in a matrix like

 (Intercept) strainA_treatmentA strainA_treatmentB  strainB_treatmentA strainB_treatmentB
       1                      0                  1                  0                   0
       1                      1                  0                  0                   0


About the error, it seems like the data is not in the right format, you have either NA or Inf somewhere in the matrix. Check it and rethink the processing you made.