Hello R users, I have a puzzle with the VGAM package, on my first excursion into generalized additive models, in that this very nice package seems to want to do either more or less than what I want. Precisely, I have a 4-component outcome, y, and am fitting multinomial logistic regression with one predictor x. What I would like to find out is, is there a single nonlinear function f(x) which acts in place of the linear predictor x. There is a mechanistic reason to believe this is sensible. So I'd like to fit a model \eta_j = \beta_{ (j) 0 } + \beta_{ (j) x } f(x) where both the function f(x) and its scaling coefficients \beta_{ (j) x } are fit simultaneously. Here \eta_j is the linear predictor, the logodds of outcome j vs the reference outcome. I cannot see how to fit exactly this. Instead I seem to be able to do the following: vgam(formula = y ~ s(x), family = multinomial) fits the model \eta_j = \beta_{ (j) 0 } + \beta_{ (j) x } f_j (x) i.e. a different function f_j (x) is fit for each outcome. vgam(formula = y ~ s(x), family = multinomial, constraints list(`(Intercept)`= diag(1,3), 's(x)' = matrix(c(1,1,1),3,1)) ) fits the model \eta_j = \beta_{ (j) 0 } + f (x) i.e. a single function f (x) is fit, but scaled the same for each outcome. I'd like one function, scaled differently. Of course, vgam(formula = y ~ x, family = multinomial) fits the model \eta_j = \beta_{ (j) 0 } + \beta_{ (j) x } x which has the scaling, but not the nonlinear function. Perhaps this is achievable using bs(), xij, and vglm, or even via the constraint matrix, but I did not succeed. Any help appreciated! Edward -- Edward Wallace, PhD Postdoctoral fellow, Drummond lab Harvard FAS center for Systems Biology ewallace at cgr.harvard.edu 773-517-4009
Edward Wallace <ewjwallace <at> gmail.com> writes:> > Hello R users, > I have a puzzle with the VGAM package, on my first excursion into > generalized additive models, in that this very nice package seems to > want to do either more or less than what I want. > > Precisely, I have a 4-component outcome, y, and am fitting multinomial > logistic regression with one predictor x. What I would like to find > out is, is there a single nonlinear function f(x) which acts in place > of the linear predictor x. There is a mechanistic reason to believe > this is sensible. So I'd like to fit a model > \eta_j = \beta_{ (j) 0 } + \beta_{ (j) x } f(x) > where both the function f(x) and its scaling coefficients \beta_{ (j) > x } are fit simultaneously. Here \eta_j is the linear predictor, the > logodds of outcome j vs the reference outcome. I cannot see how to fit > exactly this. Instead I seem to be able to do the following: >Hello, try rrvglm(y ~ 1 + bs(x), fam = multinomial, trace = TRUE) It seems what you want is a stereotype model with a smooth function. Unfortunately rrvglm() is restricted to regression splines. You could extract out the scaling coeffs and feed them into vgam() using the constraints argument, but that would not be optimal in any strict sense. cheers Thomas
Edward Wallace <ewjwallace <at> gmail.com> writes:> > Hello R users, > I have a puzzle with the VGAM package, on my first excursion into > generalized additive models, in that this very nice package seems to > want to do either more or less than what I want. > > Precisely, I have a 4-component outcome, y, and am fitting multinomial > logistic regression with one predictor x. What I would like to find > out is, is there a single nonlinear function f(x) which acts in place > of the linear predictor x. There is a mechanistic reason to believe > this is sensible. So I'd like to fit a model > \eta_j = \beta_{ (j) 0 } + \beta_{ (j) x } f(x) > where both the function f(x) and its scaling coefficients \beta_{ (j) > x } are fit simultaneously. Here \eta_j is the linear predictor, the > logodds of outcome j vs the reference outcome. I cannot see how to fit > exactly this.Thomas Yee wrote :> Hello, > > try > > rrvglm(y ~ 1 + bs(x), fam = multinomial, trace = TRUE) > > It seems what you want is a stereotype model with > a smooth function. > Unfortunately rrvglm() is restricted to regression splines.Thank you very much. This seems to work, but occasionally quits producing the cryptic error "Error in devmu[smallmu] = smy * log(smu) : NAs are not allowed in subscripted assignments" Any ideas?> You could extract out the scaling coeffs and feed them > into vgam() using the constraints argument, but that > would not be optimal in any strict sense.Really? I thought I had tried adding constraint matrices such as [1 ; 2] and vglm raised an error saying it needed entries to be 1 or 0. I can check that if you'd like. Edward -- Edward Wallace, PhD Harvard FAS center for Systems Biology +1-773-517-4009