Confidence intervals with gamlss package My question regards the use of the gamlss package. I am using gamlss package to fit a dataset to a logistic function. There is only one predictor, let me denote it with x and because the overall dependence is not exactly sigmoid, a better model is achieved by wrapping the predictor x in a smoother function. I chose the cubic spline smoother (cs). The response is binomial. I'll denote it with y:  
y = cbind(number of successful events, total number - number of successful events). 

The R code is  the following:
cs15 <- gamlss(y~cs(x, df=15), sigma.fo=~1, family=BI,  data=mydata, trace=FALSE) 

I want to estimate predicted response for a set of predictors not contained in the original data set. I know this can be achieved with the function: 
cs15fit = predict(cs25, newdata=data.frame(x=xnew), type="response")

However, my problem is that I also want to estimate the standard errors, which should be done by adding se.fit=T:
cs15fit=predict(cs25, newdata=data.frame(x=xnew), type="response", se.fit=T)

But the addition of se.fit = T produces the following error:

se.fit = TRUE is not supported for new data values at the moment

Anyone know how I can still find standard errors for the new values? 
 A: It sounds like you want prediction intervals. The se.fit argument name reflects that it will give you standard errors extracted from the original fit, and it sounds like it does this only if you do not specify newdata -- i.e. want the prediction of the original fit. (And my suspicion is that it might be giving you standard errors of your coefficients rather than prediction intervals.)
So you want to search for 'prediction interval of GAM' rather than 'standard errors'. I found a couple of papers that way, but that's just the GAM case and GAMLSS (being more flexible) is going to be even more difficult to get right.
I'd suggest falling back to a GLM -- are you sure that your "better model" isn't overfitting and that it makes a substantial difference if it's not -- where things are more straightforward and the posting @boro141 links to is applicable. (Note that you have to do more than get the se.fit results even in that case.)
http://opensiuc.lib.siu.edu/cgi/viewcontent.cgi?article=1095&context=gs_rp
A: See the answer to this question: How are the standard errors computed for the fitted values from a logistic regression?
The linked answer relates to a glm. I am not sure whether using a smoothing term as you have done affects the standard error calculation, but it should give you an idea of how standard errors are computed in predict.
