1

I have a matrix Expr with rows representing variables and columns samples. I have a categorical vector called groups (containing either "A","B", or "C") I want to test which of variables 'Expr' can be explained by the fact that the sample belong to a group.

My strategy would be modelling the problem with a generalized additive model (with a negative binomial distribution). And then I want use a likelihood ratio test in a variable wise way to get a p value for each variable. I do:

require(VGAM)
m <- vgam(Expr ~ group, family=negbinomial)
m_alternative <- vgam(Expr ~ 1, family=negbinomial)

and then:

lr <- lrtest(m, m_alternative)

The last step is wrong because it is testing the overall likelihood ratio of the two model not the variable wise. Instead of a single p value I would like to get a vector of the p-values for every variable.

How should I do it? (I am very new to R, so forgive me my stupidity)

Gioelelm
  • 2,645
  • 5
  • 30
  • 49
  • Can you explain a bit more why this isn't working? Or what you expect it to do? – Gavin Simpson May 21 '14 at 17:16
  • Basically you're testing if a model with `group` is any better than a model without it. That is an overall test for `group`. It sounds like you are doing what you want. – MrFlick May 21 '14 at 18:19
  • @MrFlick usually that is the what one wants to do but what I am trying to do here is a little bit different. I want to test for EVERY variable of Expr if the model with grouping is better than Expr ~ 1 that is basically just the noise. I will use the p-values to determine wich variable is changing between groups – Gioelelm May 21 '14 at 18:50

1 Answers1

0

It sounds like you want to use Expr as your predictors It think you may have your formula backwards. The response should be on the left, so I guess that's groups in your case.

If Expr is a data.frame, you can do regression on all variables with

m <- vgam(group ~ ., Expr, family=negbinomial)

If class(Expr)=="matrix", then 

m <- vgam(group ~ Expr, family=negbinomial)

probably should work, but you may just get slightly odd looking coefficient labels.

MrFlick
  • 195,160
  • 17
  • 277
  • 295