I fit a Generalized Additive Model in the Negative Binomial family using gam from the mgcv package. I have a data frame containing my dependent variable Y, an independent variable X, other independent variables Oth and a factor Fac. I would like to fit the following model
Y ~ s(X) + Oth
with a different theta per factor level. In other words, I use
fit <- gam(Y~s(X)+Oth, family=nb())
but this only gives me one dispersion parameter theta for the whole dataset. Instead, I believe that the mean is the same across factors, hence only one set of coefficients are required for s(X) and Oth, but the variance changes across factors, so I would like one dispersion estimate theta per level of Fac.
Naturally, fitting one model per factor level does not work because I would then get one set of coefficients for the independent variables per factor level, instead of one for the whole dataset.
The best way to solve your problem is to use a gamlss package. You will be able to model mean by X variables and variance by X and Fac variables so you will get parameters for all levels of your factor.