R : Copyright 2005, The R Foundation for Statistical Computing
Version 2.1.1  (2005-06-20), ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.

> ### * <HEADER>
> ###
> attach(NULL, name = "CheckExEnv")
> assign(".CheckExEnv", as.environment(2), pos = length(search())) # base
> ## add some hooks to label plot pages for base and grid graphics
> setHook("plot.new", ".newplot.hook")
> setHook("persp", ".newplot.hook")
> setHook("grid.newpage", ".gridplot.hook")
> 
> assign("cleanEx",
+        function(env = .GlobalEnv) {
+ 	   rm(list = ls(envir = env, all.names = TRUE), envir = env)
+            RNGkind("default", "default")
+ 	   set.seed(1)
+    	   options(warn = 1)
+ 	   delayedAssign("T", stop("T used instead of TRUE"),
+ 		  assign.env = .CheckExEnv)
+ 	   delayedAssign("F", stop("F used instead of FALSE"),
+ 		  assign.env = .CheckExEnv)
+ 	   sch <- search()
+ 	   newitems <- sch[! sch %in% .oldSearch]
+ 	   for(item in rev(newitems))
+                eval(substitute(detach(item), list(item=item)))
+ 	   missitems <- .oldSearch[! .oldSearch %in% sch]
+ 	   if(length(missitems))
+ 	       warning("items ", paste(missitems, collapse=", "),
+ 		       " have been removed from the search path")
+        },
+        env = .CheckExEnv)
> assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now
> assign("ptime", proc.time(), env = .CheckExEnv)
> grDevices::postscript("gamlss-Examples.ps")
> assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv)
> options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly"))
> options(warn = 1)    
> library('gamlss')
Loading required package: splines
 **********   GAMLSS Version 1.0-0 ********** 
For more on GAMLSS look at http://www.londonmet.ac.uk/gamlss/ 
Type gamlssNews() to see new features/changes/bug fixes.
> 
> assign(".oldSearch", search(), env = .CheckExEnv)
> assign(".oldNS", loadedNamespaces(), env = .CheckExEnv)
> cleanEx(); ..nameEx <- "BB"
> 
> ### * BB
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BB
> ### Title: Beta Binomial Distribution For Fitting a GAMLSS Model
> ### Aliases: BB dBB pBB qBB rBB
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BB()# gives information about the default links for the Beta Binomial distribution 

GAMLSS Family: BB Beta Binomial 
Link function for mu   : logit 
Link function for sigma: log 
> #plot the pdf
> plot(function(y) dBB(y, mu = .5, sigma = 1, bd =40), from=0, to=40, n=40+1, type="h")
> #calculate the cdf and plotting it
> ppBB <- pBB(seq(from=0, to=40), mu=.2 , sigma=3, bd=40)
> plot(0:40,ppBB, type="h")
> #calculating quantiles and plotting them  
> qqBB <- qBB(ppBB, mu=.2 , sigma=3, bd=40)
> plot(qqBB~ ppBB)
> # when the argument fast is useful
> p <- pBB(c(0,1,2,3,4,5), mu=.01 , sigma=1, bd=5)
> qBB(p, mu=.01 , sigma=1, bd=5, fast=TRUE)
[1] 0 1 2 3 4 5
> #  0 1 1 2 3 5
> qBB(p, mu=.01 , sigma=1, bd=5, fast=FALSE)
[1] 0 1 2 3 4 5
> #  0 1 2 3 4 5
> # generate random sample
> tN <- table(Ni <- rBB(1000, mu=.2, sigma=1, bd=20))
> r <- barplot(tN, col='lightblue')
> # fitting a model 
> data(aep)   
> library(gamlss)
> # fits a Beta-Binomial model
> h<-gamlss(y~ward+loglos+year, sigma.formula=~year+ward, family=BB, data=aep) 
GAMLSS-RS iteration 1: Global Deviance = 4490.361 
GAMLSS-RS iteration 2: Global Deviance = 4483.13 
GAMLSS-RS iteration 3: Global Deviance = 4483.021 
GAMLSS-RS iteration 4: Global Deviance = 4483.02 
GAMLSS-RS iteration 5: Global Deviance = 4483.020 
> rm(h, r, tN, ppBB, qqBB)
> 
> 
> 
> cleanEx(); ..nameEx <- "BCPE"
> 
> ### * BCPE
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BCPE
> ### Title: Box-Cox Power Exponential distribution for fitting a GAMLSS
> ### Aliases: BCPE dBCPE pBCPE qBCPE rBCPE BCPEuntr checkBCPE
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BCPE()   # 

GAMLSS Family: BCPE Box-Cox Power Exponential 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
Link function for tau  : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCPE, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4781.806 
GAMLSS-RS iteration 2: Global Deviance = 4777.83 
GAMLSS-RS iteration 3: Global Deviance = 4777.900 
GAMLSS-RS iteration 4: Global Deviance = 4777.900 
GAMLSS-RS iteration 5: Global Deviance = 4777.899 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.003533981 
                       variance   =  1.000951 
               coef. of skewness  =  0.0008490857 
               coef. of kurtosis  =  3.113359 
Filliben correlation coefficient  =  0.9989077 
*******************************************************************
> rm(h)
> plot(function(x)dBCPE(x, mu=5,sigma=.5,nu=1, tau=3), 0.0, 15, 
+  main = "The BCPE  density mu=5,sigma=.5,nu=1, tau=3")
> plot(function(x) pBCPE(x, mu=5,sigma=.5,nu=1, tau=3), 0.0, 15, 
+  main = "The BCPE  cdf mu=5, sigma=.5, nu=1, tau=3")
> 
> 
> 
> cleanEx(); ..nameEx <- "BCt"
> 
> ### * BCt
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BCT
> ### Title: Box-Cox t distribution for fitting a GAMLSS
> ### Aliases: BCT dBCT pBCT qBCT rBCT BCTuntr
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BCT()   # gives information about the default links for the Box Cox t distribution

GAMLSS Family: BCT Box-Cox t 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
Link function for tau  : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) # 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.001188135 
                       variance   =  1.002852 
               coef. of skewness  =  -0.008895612 
               coef. of kurtosis  =  2.9886 
Filliben correlation coefficient  =  0.9992773 
*******************************************************************
> plot(function(x)dBCT(x, mu=5,sigma=.5,nu=1, tau=2), 0.0, 20, 
+  main = "The BCT  density mu=5,sigma=.5,nu=1, tau=2")
> plot(function(x) pBCT(x, mu=5,sigma=.5,nu=1, tau=2), 0.0, 20, 
+  main = "The BCT  cdf mu=5, sigma=.5, nu=1, tau=2")
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "BE"
> 
> ### * BE
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BE
> ### Title: The beta distribution for fitting a GAMLSS
> ### Aliases: BE dBE pBE qBE rBE
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BE()# gives information about the default links for the normal distribution

GAMLSS Family: BE Beta 
Link function for mu   : logit 
Link function for sigma: logit 
> dat<-rBE(100, mu=.3, sigma=.5)
> hist(dat)        
> mod1<-gamlss(dat~1,family=BE) # fits a constant for mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = -58.7276 
GAMLSS-RS iteration 2: Global Deviance = -59.1803 
GAMLSS-RS iteration 3: Global Deviance = -59.1852 
GAMLSS-RS iteration 4: Global Deviance = -59.1853 
> fitted(mod1)[1]
[1] 0.3019378
> fitted(mod1,"sigma")[1]
[1] 0.4504902
> plot(function(y) dBE(y, mu=.1 ,sigma=.5), 0.001, .999)
> plot(function(y) pBE(y, mu=.1 ,sigma=.5), 0.001, 0.999)
> plot(function(y) qBE(y, mu=.1 ,sigma=.5), 0.001, 0.999)
> plot(function(y) qBE(y, mu=.1 ,sigma=.5, lower.tail=FALSE), 0.001, .999)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "BEINF"
> 
> ### * BEINF
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BEINF
> ### Title: The beta inflated distribution for fitting a GAMLSS
> ### Aliases: BEINF dBEINF pBEINF qBEINF rBEINF plotBEINF meanBEINF
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BEINF()# gives information about the default links for the normal distribution

GAMLSS Family: BEINF Beta Inflated 
Link function for mu   : logit 
Link function for sigma: logit 
Link function for nu   : log 
Link function for tau  : log 
> # plotting the distribution
> plotBEINF( mu =.5 , sigma=.5, nu = 0.5, tau = 0.5, from = 0, to=1, n = 101)
> # plotting the cdf
> plot(function(y) pBEINF(y, mu=.5 ,sigma=.5, nu = 0.5, tau = 0.5,), 0, 1)
> # plotting the inverse cdf
> plot(function(y) qBEINF(y, mu=.5 ,sigma=.5, nu = 0.5, tau = 0.5,), 0.01, .99)
> # generate random numbers
> dat <- rBEINF(1000,mu=.5,sigma=.5, nu=.5, tau=.5)
> # fit a model to the data 
> m1<-gamlss(dat~1,family=BEINF)
GAMLSS-RS iteration 1: Global Deviance = 2015.544 
GAMLSS-RS iteration 2: Global Deviance = 2012.611 
GAMLSS-RS iteration 3: Global Deviance = 2012.581 
GAMLSS-RS iteration 4: Global Deviance = 2012.581 
> meanBEINF(m1)[1]
[1] 0.499456
> 
> 
> 
> cleanEx(); ..nameEx <- "BI"
> 
> ### * BI
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BI
> ### Title: Binomial distribution for fitting a GAMLSS
> ### Aliases: BI pBI dBI qBI rBI
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BI()# gives information about the default links for the Binomial distribution 

GAMLSS Family: BI Binomial 
Link function for mu   : logit 
> data(aep)   
> library(gamlss)
> h<-gamlss(y~ward+loglos+year, family=BI, data=aep)  
GAMLSS-RS iteration 1: Global Deviance = 9452.006 
GAMLSS-RS iteration 2: Global Deviance = 9452.006 
> rm(h)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "CG"
> 
> ### * CG
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: BCCG
> ### Title: Box-Cox Cole and Green distribution (or Box-Cox normal) for
> ###   fitting a GAMLSS
> ### Aliases: BCCG BCCGuntr dBCCG pBCCG qBCCG rBCCG
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> BCCG()   # gives information about the default links for the Cole and Green distribution 

GAMLSS Family: BCCG Box-Cox-Cole-Green 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCCG, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4784.471 
GAMLSS-RS iteration 2: Global Deviance = 4783.468 
GAMLSS-RS iteration 3: Global Deviance = 4783.493 
GAMLSS-RS iteration 4: Global Deviance = 4783.494 
GAMLSS-RS iteration 5: Global Deviance = 4783.494 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.0004621194 
                       variance   =  1.001645 
               coef. of skewness  =  0.003308175 
               coef. of kurtosis  =  3.642133 
Filliben correlation coefficient  =  0.9973198 
*******************************************************************
> plot(function(x) dBCCG(x, mu=5,sigma=.5,nu=-1), 0.0, 20, 
+  main = "The BCCG  density mu=5,sigma=.5,nu=-1")
> plot(function(x) pBCCG(x, mu=5,sigma=.5,nu=-1), 0.0, 20, 
+  main = "The BCCG  cdf mu=5, sigma=.5, nu=-1")
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "GA"
> 
> ### * GA
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: GA
> ### Title: Gamma distribution for fitting a GAMLSS
> ### Aliases: GA dGA qGA pGA rGA
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> GA()# gives information about the default links for the gamma distribution      

GAMLSS Family: GA Gamma 
Link function for mu   : log 
Link function for sigma: log 
> dat<-rgamma(100, shape=1, scale=10) # generates 100 random observations 
> gamlss(dat~1,family=GA) 
GAMLSS-RS iteration 1: Global Deviance = 637.6216 
GAMLSS-RS iteration 2: Global Deviance = 637.6216 

Family:  c("GA", "Gamma") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = GA) 

Mu Coefficients:
(Intercept)  
      2.188  
Sigma Coefficients:
(Intercept)  
 -0.0008612  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     637.622 
            AIC:     641.622 
            SBC:     646.832 
> # fits a constant for each parameter mu and sigma of the gamma distribution
> newdata<-rGA(1000,mu=1,sigma=1) # generates 1000 random observations
> hist(newdata) 
> rm(dat,newdata)
> 
> 
> 
> cleanEx(); ..nameEx <- "GU"
> 
> ### * GU
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: GU
> ### Title: The Gumbel distribution for fitting a GAMLSS
> ### Aliases: GU dGU pGU qGU rGU
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> 
> plot(function(x) dGU(x, mu=0,sigma=1), -6, 3, 
+  main = "{Gumbel  density mu=0,sigma=1}")
> GU()# gives information about the default links for the Gumbel distribution      

GAMLSS Family: GU Gumbel 
Link function for mu   : identity 
Link function for sigma: log 
> dat<-rGU(100, mu=10, sigma=2) # generates 100 random observations 
> gamlss(dat~1,family=GU) # fits a constant for each parameter mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = 419.0252 
GAMLSS-RS iteration 2: Global Deviance = 419.0239 
GAMLSS-RS iteration 3: Global Deviance = 419.0239 

Family:  c("GU", "Gumbel") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = GU) 

Mu Coefficients:
(Intercept)  
      10.05  
Sigma Coefficients:
(Intercept)  
     0.5196  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     419.024 
            AIC:     423.024 
            SBC:     428.234 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "IC"
> 
> ### * IC
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: IC
> ### Title: Gives the GAIC for a GAMLSS Object
> ### Aliases: IC AIC.gamlss GAIC extractAIC.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> mod1<-gamlss(y~cs(x,df=3),sigma.fo=~cs(x,df=1),family=BCT, data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> IC(mod1)
[1] 4795.878
> mod2<-gamlss(y~cs(x,df=3),sigma.fo=~x,family=BCT, data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4789.114 
GAMLSS-RS iteration 2: Global Deviance = 4788.652 
GAMLSS-RS iteration 3: Global Deviance = 4788.64 
GAMLSS-RS iteration 4: Global Deviance = 4788.639 
GAMLSS-RS iteration 5: Global Deviance = 4788.639 
> AIC(mod1,mod2,k=3)
            df      AIC
mod1 10.000899 4805.879
mod2  8.999233 4815.636
> GAIC(mod1,mod2,k=3)
            df      AIC
mod1 10.000899 4805.879
mod2  8.999233 4815.636
> extractAIC(mod1,k=3)
[1]   10.0009 4805.8786
> rm(mod1,mod2)
> 
> 
> 
> cleanEx(); ..nameEx <- "IG"
> 
> ### * IG
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: IG
> ### Title: Inverse Gaussian distribution for fitting a GAMLSS
> ### Aliases: IG dIG pIG qIG rIG
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> IG()# gives information about the default links for the normal distribution

GAMLSS Family: IG Inverse Gaussian 
Link function for mu   : log 
Link function for sigma: log 
> data(rent)        
> gamlss(R~cs(Fl),family=IG, data=rent) # 
GAMLSS-RS iteration 1: Global Deviance = 28268.68 
GAMLSS-RS iteration 2: Global Deviance = 28268.68 

Family:  c("IG", "Inverse Gaussian") 
Fitting method: RS() 

Call:  gamlss(formula = R ~ cs(Fl), family = IG, data = rent) 

Mu Coefficients:
(Intercept)       cs(Fl)  
    5.95648      0.01062  
Sigma Coefficients:
(Intercept)  
     -4.126  

Degrees of Freedom for the fit: 6.000595 Residual Deg. of Freedom   1962.999 
Global Deviance:     28268.7 
            AIC:     28280.7 
            SBC:     28314.2 
> plot(function(x)dIG(x, mu=1,sigma=.5), 0.01, 6, 
+  main = "{Inverse Gaussian  density mu=1,sigma=0.5}")
> plot(function(x)pIG(x, mu=1,sigma=.5), 0.01, 6, 
+  main = "{Inverse Gaussian  cdf mu=1,sigma=0.5}")
> 
> 
> 
> cleanEx(); ..nameEx <- "JSU"
> 
> ### * JSU
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: JSU
> ### Title: The Johnson's Su distribution for fitting a GAMLSS
> ### Aliases: JSU dJSU pJSU qJSU rJSU
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> JSU()   

GAMLSS Family: JSU Johnson SU 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
Link function for tau  : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=JSU, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4775.214 
GAMLSS-RS iteration 2: Global Deviance = 4775.664 
GAMLSS-RS iteration 3: Global Deviance = 4775.803 
GAMLSS-RS iteration 4: Global Deviance = 4775.841 
GAMLSS-RS iteration 5: Global Deviance = 4775.852 
GAMLSS-RS iteration 6: Global Deviance = 4775.856 
GAMLSS-RS iteration 7: Global Deviance = 4775.857 
GAMLSS-RS iteration 8: Global Deviance = 4775.858 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  -0.0009393478 
                       variance   =  1.001566 
               coef. of skewness  =  -0.01301643 
               coef. of kurtosis  =  2.986125 
Filliben correlation coefficient  =  0.99931 
*******************************************************************
> rm(h)
> plot(function(x)dJSU(x, mu=0,sigma=1,nu=-1, tau=.5), -4, 4, 
+  main = "The JSU  density mu=0,sigma=1,nu=-1, tau=.5")
> plot(function(x) pJSU(x, mu=0,sigma=1,nu=-1, tau=.5), -4, 4, 
+  main = "The JSU  cdf mu=0, sigma=1, nu=-1, tau=.5")
> 
> 
> 
> cleanEx(); ..nameEx <- "JSUoriginal"
> 
> ### * JSUoriginal
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: JSUo
> ### Title: The original Johnson's Su distribution for fitting a GAMLSS
> ### Aliases: JSUo dJSUo pJSUo qJSUo rJSUo
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> JSU()   

GAMLSS Family: JSU Johnson SU 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
Link function for tau  : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=JSUo, 
+           data=abdom, method=mixed(2,20)) 
GAMLSS-RS iteration 1: Global Deviance = 5399.72 
GAMLSS-RS iteration 2: Global Deviance = 4910.275 
GAMLSS-CG iteration 1: Global Deviance = 4783.778 
GAMLSS-CG iteration 2: Global Deviance = 4776.15 
GAMLSS-CG iteration 3: Global Deviance = 4775.891 
GAMLSS-CG iteration 4: Global Deviance = 4775.821 
GAMLSS-CG iteration 5: Global Deviance = 4775.822 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.0001633889 
                       variance   =  1.000783 
               coef. of skewness  =  0.01226295 
               coef. of kurtosis  =  2.968972 
Filliben correlation coefficient  =  0.9993104 
*******************************************************************
> rm(h)
> plot(function(x)dJSUo(x, mu=0,sigma=1,nu=-1, tau=.5), -4, 15, 
+  main = "The JSUo  density mu=0,sigma=1,nu=-1, tau=.5")
> plot(function(x) pJSUo(x, mu=0,sigma=1,nu=-1, tau=.5), -4, 15, 
+  main = "The JSUo  cdf mu=0, sigma=1, nu=-1, tau=.5")
> 
> 
> 
> cleanEx(); ..nameEx <- "LNO"
> 
> ### * LNO
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: LNO
> ### Title: Log Normal distribution for fitting a GAMLSS
> ### Aliases: LNO dLNO pLNO qLNO rLNO
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> LNO()# gives information about the default links for the Box Cox distribution 

GAMLSS Family: LNO Log Normal 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : 
> data(abdom)
> h1<-gamlss(y~cs(x), family=LNO, data=abdom)#fits the log-Normal distribution   
GAMLSS-RS iteration 1: Global Deviance = 4870.896 
GAMLSS-RS iteration 2: Global Deviance = 4870.896 
> # to change to square root transformation, i.e. fix nu=0.5 
> h2<-gamlss(y~cs(x), family=LNO, data=abdom, nu.fix=TRUE, nu.start=0.5)
GAMLSS-RS iteration 1: Global Deviance = 4826.368 
GAMLSS-RS iteration 2: Global Deviance = 4826.368 
> rm(h1,h2)
> 
> 
> 
> cleanEx(); ..nameEx <- "LOGISTIC"
> 
> ### * LOGISTIC
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: LO
> ### Title: Logistic distribution for fitting a GAMLSS
> ### Aliases: LO dLO pLO qLO rLO
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> LO()# gives information about the default links for the Logistic distribution 

GAMLSS Family: LO Logistic 
Link function for mu   : identity 
Link function for sigma: log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=LO, data=abdom) # fits 
GAMLSS-RS iteration 1: Global Deviance = 4781.45 
GAMLSS-RS iteration 2: Global Deviance = 4781.096 
GAMLSS-RS iteration 3: Global Deviance = 4781.118 
GAMLSS-RS iteration 4: Global Deviance = 4781.118 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.01518417 
                       variance   =  0.9897921 
               coef. of skewness  =  0.1511567 
               coef. of kurtosis  =  2.784563 
Filliben correlation coefficient  =  0.9982932 
*******************************************************************
> rm(h)
> plot(function(y) dLO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) pLO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) qLO(y, mu=10 ,sigma=2), 0, 1)
> 
> 
> 
> cleanEx(); ..nameEx <- "Mums"
> 
> ### * Mums
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Mums
> ### Title: Mothers encouragement data
> ### Aliases: Mums
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(Mums)
> MM<-xtabs(~mums+qual, data=Mums)
> mosaicplot(MM, color=TRUE)
> MM<-xtabs(~mums+ethn+gender, data=Mums)
> mosaicplot(MM, color=TRUE)
> 
> 
> 
> cleanEx(); ..nameEx <- "NBI"
> 
> ### * NBI
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: NBI
> ### Title: Negative Binomial type I distribution for fitting a GAMLSS
> ### Aliases: NBI dNBI pNBI qNBI rNBI
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> NBI()   # gives information about the default links for the Negative Binomial type I distribution  

GAMLSS Family: NBI Negative Binomial type I 
Link function for mu   : log 
Link function for sigma: log 
> data(aids)
> h<-gamlss(y~cs(x,df=7)+qrt, family=NBI, data=aids) # fits the model 
GAMLSS-RS iteration 1: Global Deviance = 365.8195 
GAMLSS-RS iteration 2: Global Deviance = 361.9924 
GAMLSS-RS iteration 3: Global Deviance = 362.1084 
GAMLSS-RS iteration 4: Global Deviance = 362.1114 
GAMLSS-RS iteration 5: Global Deviance = 362.1123 
> plot(h)
*******************************************************************
	 Summary of the Randomised Quantile Residuals
                           mean   =  -0.01086206 
                       variance   =  0.9525586 
               coef. of skewness  =  -0.5769322 
               coef. of kurtosis  =  3.25145 
Filliben correlation coefficient  =  0.9853362 
*******************************************************************
> # plotting the distribution
> plot(function(y) dNBI(y, mu = 10, sigma = 0.5 ), from=0, to=40, n=40+1, type="h")
> pdf.plot(family=NBI, mu=10, sigma=0.5, min=0, max=40, step=1)
> # creating random variables and plot them 
> tN <- table(Ni <- rNBI(1000, mu=5, sigma=0.5))
> r <- barplot(tN, col='lightblue')
> rm(h, tN, r)
> 
> 
> 
> cleanEx(); ..nameEx <- "NBII"
> 
> ### * NBII
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: NBII
> ### Title: Negative Binomial type II distribution for fitting a GAMLSS
> ### Aliases: NBII dNBII pNBII qNBII rNBII
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> NBII()  # gives information about the default links for the Negative Binomial type II distribution  

GAMLSS Family: NBII Negative Binomial type II 
Link function for mu   : log 
Link function for sigma: log 
> data(aids)
> h<-gamlss(y~cs(x,df=7)+qrt, family=NBII, data=aids) # fits a model 
GAMLSS-RS iteration 1: Global Deviance = 466.8198 
GAMLSS-RS iteration 2: Global Deviance = 395.6372 
GAMLSS-RS iteration 3: Global Deviance = 366.6319 
GAMLSS-RS iteration 4: Global Deviance = 365.9296 
GAMLSS-RS iteration 5: Global Deviance = 365.9567 
GAMLSS-RS iteration 6: Global Deviance = 365.9569 
> plot(h)
*******************************************************************
	 Summary of the Randomised Quantile Residuals
                           mean   =  -0.0162118 
                       variance   =  1.025826 
               coef. of skewness  =  -0.5469316 
               coef. of kurtosis  =  3.761933 
Filliben correlation coefficient  =  0.9830424 
*******************************************************************
> # plotting the distribution
> plot(function(y) dNBII(y, mu = 10, sigma = 0.5 ), from=0, to=40, n=40+1, type="h")
> pdf.plot(family=NBII, mu=10, sigma=0.5, min=0, max=40, step=1)
> # creating random variables and plot them 
> tN <- table(Ni <- rNBII(1000, mu=5, sigma=0.5))
> r <- barplot(tN, col='lightblue')
> rm(h, tN, r)
> 
> rm(h)
Warning: remove: variable "h" was not found
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "NET"
> 
> ### * NET
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: NET
> ### Title: Normal Exponential t distribution (NET) for fitting a GAMLSS
> ### Aliases: NET dNET pNET
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> NET()   # 

GAMLSS Family: NET Normal Exponential t 
Link function for mu   : 
Link function for sigma: 
Link function for nu   : 
Link function for tau  : 
> data(abdom)
> # fit NET with nu=1 and tau=3
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=NET, 
+         data=abdom, nu.start=2, tau.start=3) 
GAMLSS-RS iteration 1: Global Deviance = 4780.713 
GAMLSS-RS iteration 2: Global Deviance = 4779.658 
GAMLSS-RS iteration 3: Global Deviance = 4779.681 
GAMLSS-RS iteration 4: Global Deviance = 4779.681 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.0002559666 
                       variance   =  0.9886621 
               coef. of skewness  =  0.1127455 
               coef. of kurtosis  =  2.87434 
Filliben correlation coefficient  =  0.9983134 
*******************************************************************
> rm(h)
> plot(function(x)dNET(x, mu=0,sigma=1,nu=2, tau=3), -5, 5)
> plot(function(x)pNET(x, mu=0,sigma=1,nu=2, tau=3), -5, 5) 
> 
> 
> 
> cleanEx(); ..nameEx <- "NO"
> 
> ### * NO
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: NO
> ### Title: Normal distribution for fitting a GAMLSS
> ### Aliases: NO dNO pNO qNO rNO
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> NO()# gives information about the default links for the normal distribution

GAMLSS Family: NO Normal 
Link function for mu   : identity 
Link function for sigma: log 
> dat<-rNO(100)        
> gamlss(dat~1,family=NO) # fits a constant for mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = 261.31 
GAMLSS-RS iteration 2: Global Deviance = 261.31 

Family:  c("NO", "Normal") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = NO) 

Mu Coefficients:
(Intercept)  
     0.1089  
Sigma Coefficients:
(Intercept)  
    -0.1124  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     261.31 
            AIC:     265.31 
            SBC:     270.52 
> plot(function(y) dNO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) pNO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) qNO(y, mu=10 ,sigma=2), 0, 1)
> 
> 
> 
> cleanEx(); ..nameEx <- "NOvar"
> 
> ### * NOvar
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: NO.var
> ### Title: Normal distribution (with variance as sigma parameter) for
> ###   fitting a GAMLSS
> ### Aliases: NO.var dNO.var pNO.var qNO.var rNO.var
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> NO()# gives information about the default links for the normal distribution

GAMLSS Family: NO Normal 
Link function for mu   : identity 
Link function for sigma: log 
> dat<-rNO(100)        
> gamlss(dat~1,family=NO) # fits a constant for mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = 261.31 
GAMLSS-RS iteration 2: Global Deviance = 261.31 

Family:  c("NO", "Normal") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = NO) 

Mu Coefficients:
(Intercept)  
     0.1089  
Sigma Coefficients:
(Intercept)  
    -0.1124  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     261.31 
            AIC:     265.31 
            SBC:     270.52 
> plot(function(y) dNO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) pNO(y, mu=10 ,sigma=2), 0, 20)
> plot(function(y) qNO(y, mu=10 ,sigma=2), 0, 1)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "PE"
> 
> ### * PE
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: PE
> ### Title: Power Exponential distribution for fitting a GAMLSS
> ### Aliases: PE dPE pPE qPE rPE
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> PE()# gives information about the default links for the Power Exponential distribution  

GAMLSS Family: PE Power Exponential 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=PE, data=abdom) # fits 
GAMLSS-RS iteration 1: Global Deviance = 4794.24 
GAMLSS-RS iteration 2: Global Deviance = 4781.199 
GAMLSS-RS iteration 3: Global Deviance = 4781.230 
GAMLSS-RS iteration 4: Global Deviance = 4781.231 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.01031142 
                       variance   =  1.000830 
               coef. of skewness  =  0.1874644 
               coef. of kurtosis  =  3.117730 
Filliben correlation coefficient  =  0.998056 
*******************************************************************
> rm(h)
> # leptokurtotic
> plot(function(x) dPE(x, mu=10,sigma=2,nu=1), 0.0, 20, 
+  main = "The PE  density mu=10,sigma=2,nu=1")
> # platykurtotic
> plot(function(x) dPE(x, mu=10,sigma=2,nu=4), 0.0, 20, 
+  main = "The PE  density mu=10,sigma=2,nu=4") 
> 
> 
> 
> cleanEx(); ..nameEx <- "PIG"
> 
> ### * PIG
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: PIG
> ### Title: The Poisson-inverse Gaussian distribution for fitting a GAMLSS
> ###   model
> ### Aliases: PIG dPIG pPIG qPIG rPIG
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> PIG()# gives information about the default links for the  Poisson-inverse Gaussian distribution 

GAMLSS Family: PIG Poisson.Inverse.Gaussian 
Link function for mu   : log 
Link function for sigma: log 
> #plot the pdf using plot 
> plot(function(y) dPIG(y, mu=10, sigma = 1 ), from=0, to=50, n=50+1, type="h") # pdf
> # plot the cdf
> plot(seq(from=0,to=50),pPIG(seq(from=0,to=50), mu=10, sigma=1), type="h")   # cdf
> # generate random sample
> tN <- table(Ni <- rPIG(1000, mu=5, sigma=1, fast=TRUE))
> r <- barplot(tN, col='lightblue')
> # fit a model to the data 
> gamlss(Ni~1,family=PIG)
GAMLSS-RS iteration 1: Global Deviance = 5069.145 
GAMLSS-RS iteration 2: Global Deviance = 5067.969 
GAMLSS-RS iteration 3: Global Deviance = 5067.92 
GAMLSS-RS iteration 4: Global Deviance = 5067.917 
GAMLSS-RS iteration 5: Global Deviance = 5067.916 
GAMLSS-RS iteration 6: Global Deviance = 5067.917 

Family:  c("PIG", "Poisson.Inverse.Gaussian") 
Fitting method: RS() 

Call:  gamlss(formula = Ni ~ 1, family = PIG) 

Mu Coefficients:
(Intercept)  
      1.430  
Sigma Coefficients:
(Intercept)  
     0.7435  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   998 
Global Deviance:     5067.92 
            AIC:     5071.92 
            SBC:     5081.73 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "PO"
> 
> ### * PO
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: PO
> ### Title: Poisson distribution for fitting a GAMLSS model
> ### Aliases: PO dPO pPO qPO rPO
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> PO()# gives information about the default links for the Poisson distribution  

GAMLSS Family: PO Poisson 
Link function for mu   : log 
> # fitting data using PO()
> data(aids)
> h<-gamlss(y~cs(x,df=7)+qrt, family=PO, data=aids) # fits the constant+x+qrt model 
GAMLSS-RS iteration 1: Global Deviance = 374.1748 
GAMLSS-RS iteration 2: Global Deviance = 374.1748 
> plot(h)
*******************************************************************
	 Summary of the Randomised Quantile Residuals
                           mean   =  -0.02949705 
                       variance   =  1.776076 
               coef. of skewness  =  -0.5494088 
               coef. of kurtosis  =  3.766308 
Filliben correlation coefficient  =  0.9831796 
*******************************************************************
> # plotting the distribution
> plot(function(y) dPO(y, mu=10 ), from=0, to=20, n=20+1, type="h")
> pdf.plot(family=PO, mu=10, min=0, max=20, step=1)
> # creating random variables and plot them 
> tN <- table(Ni <- rPO(1000, mu=5))
>  r <- barplot(tN, col='lightblue')
> rm(h, tN, r)
> 
> 
> 
> cleanEx(); ..nameEx <- "Q.stats"
> 
> ### * Q.stats
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Q.stats
> ### Title: A function to calculate the Q-statistics
> ### Aliases: Q.stats
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> Q.stats(h,xvar=abdom$x,n.inter=8)
                          Z1          Z2           Z3          Z4 AgostinoK2
12.22 to 16.78     0.2418782  0.76665696 -0.113691950 -0.83804044  0.7152376
16.78 to 20.07     0.5541892 -0.32506259 -0.359052265  0.28270686  0.2088417
20.07 to 23.5     -0.3758341 -0.18980830 -0.004821328 -2.38144971  5.6713260
23.5 to 27.07     -0.5784963 -1.45332910  1.745247716  1.83285144  6.4052340
27.07 to 30.64    -0.4425575  0.01361803  0.369817944 -0.55935775  0.4496464
30.64 to 34.5      1.1491899  0.68238279 -0.186924502  0.84860204  0.7550662
34.5 to 38.36     -0.1112958  0.42791425 -0.599987531  0.83618078  1.0591833
38.36 to 42.5     -0.4053002  0.15956057 -0.864809981  1.01942877  1.7871313
TOTAL Q stats      2.5346896  3.51602317  4.467344647 12.58432191 17.0516666
df for Q stats     2.9991766  5.99996204  7.000000000  7.00000000 14.0000000
p-val for Q stats  0.4689038  0.74183264  0.724645880  0.08290769  0.2534350
                    N
12.22 to 16.78     79
16.78 to 20.07     75
20.07 to 23.5      76
23.5 to 27.07      77
27.07 to 30.64     74
30.64 to 34.5      82
34.5 to 38.36      71
38.36 to 42.5      76
TOTAL Q stats     610
df for Q stats      0
p-val for Q stats   0
> Q.stats(h,xvar=abdom$x,n.inter=8,zvals=FALSE)
                          Q1           Q2           Q3          Q4 AgostinoK2
12.22 to 16.78    0.05850507 0.5877628961 1.292586e-02  0.70231178  0.7152376
16.78 to 20.07    0.30712571 0.1056656870 1.289185e-01  0.07992317  0.2088417
20.07 to 23.5     0.14125129 0.0360271926 2.324520e-05  5.67130273  5.6713260
23.5 to 27.07     0.33465800 2.1121654871 3.045890e+00  3.35934440  6.4052340
27.07 to 30.64    0.19585714 0.0001854507 1.367653e-01  0.31288110  0.4496464
30.64 to 34.5     1.32063743 0.4656462771 3.494077e-02  0.72012542  0.7550662
34.5 to 38.36     0.01238676 0.1831106024 3.599850e-01  0.69919829  1.0591833
38.36 to 42.5     0.16426824 0.0254595752 7.478963e-01  1.03923501  1.7871313
TOTAL Q stats     2.53468965 3.5160231683 4.467345e+00 12.58432191 17.0516666
df for Q stats    2.99917659 5.9999620431 7.000000e+00  7.00000000 14.0000000
p-val for Q stats 0.46890379 0.7418326383 7.246459e-01  0.08290769  0.2534350
                    N
12.22 to 16.78     79
16.78 to 20.07     75
20.07 to 23.5      76
23.5 to 27.07      77
27.07 to 30.64     74
30.64 to 34.5      82
34.5 to 38.36      71
38.36 to 42.5      76
TOTAL Q stats     610
df for Q stats      0
p-val for Q stats   0
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "RG"
> 
> ### * RG
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: RG
> ### Title: The Reverse Gumbel distribution for fitting a GAMLSS
> ### Aliases: RG dRG pRG qRG rRG
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> 
> plot(function(x) dRG(x, mu=0,sigma=1), -3, 6, 
+  main = "{Reverse Gumbel  density mu=0,sigma=1}")
> RG()# gives information about the default links for the Gumbel distribution      

GAMLSS Family: RG Reverse Gumbel 
Link function for mu   : identity 
Link function for sigma: log 
> dat<-rRG(100, mu=10, sigma=2) # generates 100 random observations 
> gamlss(dat~1,family=RG) # fits a constant for each parameter mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = 427.1582 
GAMLSS-RS iteration 2: Global Deviance = 427.1059 
GAMLSS-RS iteration 3: Global Deviance = 427.1054 

Family:  c("RG", "Reverse Gumbel") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = RG) 

Mu Coefficients:
(Intercept)  
      10.17  
Sigma Coefficients:
(Intercept)  
     0.5738  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     427.105 
            AIC:     431.105 
            SBC:     436.316 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "SEP"
> 
> ### * SEP
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: SEP
> ### Title: The Skew Power exponential (SEP) distribution for fitting a
> ###   GAMLSS
> ### Aliases: SEP dSEP pSEP qSEP rSEP
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> SEP()   # 

GAMLSS Family: SEP Skew Exponential Power 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : identity 
Link function for tau  : log 
> plot(function(x)dSEP(x, mu=0,sigma=1, nu=1, tau=2), -5, 5, 
+  main = "The SEP  density mu=0,sigma=1,nu=1, tau=2")
> plot(function(x) pSEP(x, mu=0,sigma=1,nu=1, tau=2), -5, 5, 
+  main = "The BCPE  cdf mu=0, sigma=1, nu=1, tau=2")
> dat<-rSEP(100,mu=10,sigma=1,nu=-1,tau=1.5)
> gamlss(dat~1,family=SEP, control=gamlss.control(n.cyc=30))
GAMLSS-RS iteration 1: Global Deviance = 224.8745 
GAMLSS-RS iteration 2: Global Deviance = 224.8387 
GAMLSS-RS iteration 3: Global Deviance = 224.8132 
GAMLSS-RS iteration 4: Global Deviance = 224.7903 
GAMLSS-RS iteration 5: Global Deviance = 224.768 
GAMLSS-RS iteration 6: Global Deviance = 224.7467 
GAMLSS-RS iteration 7: Global Deviance = 224.7255 
GAMLSS-RS iteration 8: Global Deviance = 224.704 
GAMLSS-RS iteration 9: Global Deviance = 224.6816 
GAMLSS-RS iteration 10: Global Deviance = 224.6579 
GAMLSS-RS iteration 11: Global Deviance = 224.6332 
GAMLSS-RS iteration 12: Global Deviance = 224.609 
GAMLSS-RS iteration 13: Global Deviance = 224.5826 
GAMLSS-RS iteration 14: Global Deviance = 224.5568 
GAMLSS-RS iteration 15: Global Deviance = 224.532 
GAMLSS-RS iteration 16: Global Deviance = 224.5062 
GAMLSS-RS iteration 17: Global Deviance = 224.4784 
GAMLSS-RS iteration 18: Global Deviance = 224.4463 
GAMLSS-RS iteration 19: Global Deviance = 224.4105 
GAMLSS-RS iteration 20: Global Deviance = 224.3761 
GAMLSS-RS iteration 21: Global Deviance = 224.3386 
GAMLSS-RS iteration 22: Global Deviance = 224.2925 
GAMLSS-RS iteration 23: Global Deviance = 224.2425 
GAMLSS-RS iteration 24: Global Deviance = 224.1922 
GAMLSS-RS iteration 25: Global Deviance = 224.1338 
GAMLSS-RS iteration 26: Global Deviance = 224.0522 
GAMLSS-RS iteration 27: Global Deviance = 223.9427 
GAMLSS-RS iteration 28: Global Deviance = 223.8395 
GAMLSS-RS iteration 29: Global Deviance = 223.7065 
GAMLSS-RS iteration 30: Global Deviance = 223.512 
Warning in RS() : Algorithm RS has not yet converged

Family:  c("SEP", "Skew Exponential Power") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = SEP, control = gamlss.control(n.cyc = 30)) 

Mu Coefficients:
(Intercept)  
      9.683  
Sigma Coefficients:
(Intercept)  
    -0.3168  
Nu Coefficients:
(Intercept)  
    -0.4656  
Tau Coefficients:
(Intercept)  
     0.4673  

Degrees of Freedom for the fit: 4 Residual Deg. of Freedom   96 
Global Deviance:     223.512 
            AIC:     231.512 
            SBC:     241.933 
> 
> 
> 
> cleanEx(); ..nameEx <- "SI"
> 
> ### * SI
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: SI
> ### Title: The Sichel dustribution for fitting a GAMLSS model
> ### Aliases: SI dSI pSI qSI rSI
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> SI()# gives information about the default links for the  Sichel distribution 

GAMLSS Family: SI Sichel 
Link function for mu   : log 
Link function for sigma: log 
Link function for nu   : identity 
> #plot the pdf using plot 
> plot(function(y) dSI(y, mu=10, sigma=1, nu=1), from=0, to=100, n=100+1, type="h") # pdf
> # plot the cdf
> plot(seq(from=0,to=100),pSI(seq(from=0,to=100), mu=10, sigma=1, nu=1), type="h")   # cdf
> # generate random sample
> tN <- table(Ni <- rSI(1000, mu=5, sigma=1, nu=1, fast=TRUE))
> r <- barplot(tN, col='lightblue')
> # fit a model to the data 
> gamlss(Ni~1,family=SI, control=gamlss.control(n.cyc=50))
GAMLSS-RS iteration 1: Global Deviance = 7148.412 
GAMLSS-RS iteration 2: Global Deviance = 7146.443 
GAMLSS-RS iteration 3: Global Deviance = 7145.184 
GAMLSS-RS iteration 4: Global Deviance = 7144.009 
GAMLSS-RS iteration 5: Global Deviance = 7142.855 
GAMLSS-RS iteration 6: Global Deviance = 7141.746 
GAMLSS-RS iteration 7: Global Deviance = 7140.652 
GAMLSS-RS iteration 8: Global Deviance = 7139.604 
GAMLSS-RS iteration 9: Global Deviance = 7138.59 
GAMLSS-RS iteration 10: Global Deviance = 7137.594 
GAMLSS-RS iteration 11: Global Deviance = 7136.636 
GAMLSS-RS iteration 12: Global Deviance = 7135.697 
GAMLSS-RS iteration 13: Global Deviance = 7134.8 
GAMLSS-RS iteration 14: Global Deviance = 7133.928 
GAMLSS-RS iteration 15: Global Deviance = 7133.093 
GAMLSS-RS iteration 16: Global Deviance = 7132.277 
GAMLSS-RS iteration 17: Global Deviance = 7131.487 
GAMLSS-RS iteration 18: Global Deviance = 7130.734 
GAMLSS-RS iteration 19: Global Deviance = 7129.994 
GAMLSS-RS iteration 20: Global Deviance = 7129.288 
GAMLSS-RS iteration 21: Global Deviance = 7128.598 
GAMLSS-RS iteration 22: Global Deviance = 7127.922 
GAMLSS-RS iteration 23: Global Deviance = 7127.272 
GAMLSS-RS iteration 24: Global Deviance = 7126.64 
GAMLSS-RS iteration 25: Global Deviance = 7126.038 
GAMLSS-RS iteration 26: Global Deviance = 7125.47 
GAMLSS-RS iteration 27: Global Deviance = 7124.901 
GAMLSS-RS iteration 28: Global Deviance = 7124.35 
GAMLSS-RS iteration 29: Global Deviance = 7123.829 
GAMLSS-RS iteration 30: Global Deviance = 7123.314 
GAMLSS-RS iteration 31: Global Deviance = 7122.824 
GAMLSS-RS iteration 32: Global Deviance = 7122.357 
GAMLSS-RS iteration 33: Global Deviance = 7121.894 
GAMLSS-RS iteration 34: Global Deviance = 7121.454 
GAMLSS-RS iteration 35: Global Deviance = 7121.017 
GAMLSS-RS iteration 36: Global Deviance = 7120.596 
GAMLSS-RS iteration 37: Global Deviance = 7120.192 
GAMLSS-RS iteration 38: Global Deviance = 7119.808 
GAMLSS-RS iteration 39: Global Deviance = 7119.428 
GAMLSS-RS iteration 40: Global Deviance = 7119.085 
GAMLSS-RS iteration 41: Global Deviance = 7118.73 
GAMLSS-RS iteration 42: Global Deviance = 7118.392 
GAMLSS-RS iteration 43: Global Deviance = 7118.06 
GAMLSS-RS iteration 44: Global Deviance = 7117.742 
GAMLSS-RS iteration 45: Global Deviance = 7117.445 
GAMLSS-RS iteration 46: Global Deviance = 7117.149 
GAMLSS-RS iteration 47: Global Deviance = 7116.873 
GAMLSS-RS iteration 48: Global Deviance = 7116.594 
GAMLSS-RS iteration 49: Global Deviance = 7116.325 
GAMLSS-RS iteration 50: Global Deviance = 7116.082 
Warning in RS() : Algorithm RS has not yet converged

Family:  c("SI", "Sichel") 
Fitting method: RS() 

Call:  gamlss(formula = Ni ~ 1, family = SI, control = gamlss.control(n.cyc = 50)) 

Mu Coefficients:
(Intercept)  
      1.741  
Sigma Coefficients:
(Intercept)  
     0.2071  
Nu Coefficients:
(Intercept)  
     0.5026  

Degrees of Freedom for the fit: 3 Residual Deg. of Freedom   997 
Global Deviance:     7116.08 
            AIC:     7122.08 
            SBC:     7136.81 
> 
> 
> 
> cleanEx(); ..nameEx <- "TF"
> 
> ### * TF
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: TF
> ### Title: t family distribution for fitting a GAMLSS
> ### Aliases: TF dTF pTF qTF rTF
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> TF()# gives information about the default links for the t-family distribution 

GAMLSS Family: TF t Family 
Link function for mu   : identity 
Link function for sigma: log 
Link function for nu   : log 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=TF, data=abdom) # fits 
GAMLSS-RS iteration 1: Global Deviance = 4779.744 
GAMLSS-RS iteration 2: Global Deviance = 4779.158 
GAMLSS-RS iteration 3: Global Deviance = 4779.165 
GAMLSS-RS iteration 4: Global Deviance = 4779.164 
GAMLSS-RS iteration 5: Global Deviance = 4779.163 
> plot(h)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.009694789 
                       variance   =  1.002052 
               coef. of skewness  =  0.1678607 
               coef. of kurtosis  =  2.985408 
Filliben correlation coefficient  =  0.9984764 
*******************************************************************
> newdata<-rTF(1000,mu=0,sigma=1,nu=5) # generates 1000 random observations
> hist(newdata) 
> rm(h,newdata)
> 
> 
> 
> cleanEx(); ..nameEx <- "WEI"
> 
> ### * WEI
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: WEI
> ### Title: Weibull distribution for fitting a GAMLSS
> ### Aliases: WEI dWEI pWEI qWEI rWEI
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> WEI()

GAMLSS Family: WEI Weibull 
Link function for mu   : log 
Link function for sigma: log 
> dat<-rWEI(100, mu=10, sigma=2)
> gamlss(dat~1, family=WEI)
GAMLSS-RS iteration 1: Global Deviance = 552.117 
GAMLSS-RS iteration 2: Global Deviance = 552.0641 
GAMLSS-RS iteration 3: Global Deviance = 552.0636 

Family:  c("WEI", "Weibull") 
Fitting method: RS() 

Call:  gamlss(formula = dat ~ 1, family = WEI) 

Mu Coefficients:
(Intercept)  
      2.260  
Sigma Coefficients:
(Intercept)  
     0.8125  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     552.064 
            AIC:     556.064 
            SBC:     561.274 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "WEI2"
> 
> ### * WEI2
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: WEI2
> ### Title: A specific parameterization of the Weibull distribution for
> ###   fitting a GAMLSS
> ### Aliases: WEI2 dWEI2 pWEI2 qWEI2 rWEI2
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> WEI2()

GAMLSS Family: WEI2 Weibull 
Link function for mu   : log 
Link function for sigma: log 
> dat<-rWEI(100, mu=.1, sigma=2)
> gamlss(dat~1, family=WEI2, method=CG())
GAMLSS-CG iteration 1: Global Deviance = -352.0256 
GAMLSS-CG iteration 2: Global Deviance = -367.5879 
GAMLSS-CG iteration 3: Global Deviance = -368.9265 
GAMLSS-CG iteration 4: Global Deviance = -368.9641 
GAMLSS-CG iteration 5: Global Deviance = -368.965 

Family:  c("WEI2", "Weibull") 
Fitting method: CG() 

Call:  gamlss(formula = dat ~ 1, family = WEI2, method = CG()) 

Mu Coefficients:
(Intercept)  
      5.255  
Sigma Coefficients:
(Intercept)  
      0.807  

Degrees of Freedom for the fit: 2 Residual Deg. of Freedom   98 
Global Deviance:     -368.965 
            AIC:     -364.965 
            SBC:     -359.755 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "ZIP"
> 
> ### * ZIP
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ZIP
> ### Title: Zero inflated poisson distribution for fitting a GAMLSS model
> ### Aliases: ZIP dZIP pZIP qZIP rZIP
> ### Keywords: distribution regression
> 
> ### ** Examples
> 
> ZIP()# gives information about the default links for the normal distribution

GAMLSS Family: PZI Poisson Zero Inflated 
Link function for mu   : log 
Link function for sigma: logit 
> # creating data and plotting them 
> dat<-rZIP(1000, mu=5, sigma=.1)
> r <- barplot(table(dat), col='lightblue')
> # fit the disteibution 
> mod1<-gamlss(dat~1, family=ZIP)# fits a constant for mu and sigma 
GAMLSS-RS iteration 1: Global Deviance = 4592.071 
GAMLSS-RS iteration 2: Global Deviance = 4591.955 
GAMLSS-RS iteration 3: Global Deviance = 4591.955 
> fitted(mod1)[1]
[1] 4.982282
> fitted(mod1,"sigma")[1]
[1] 0.09479173
> 
> 
> 
> cleanEx(); ..nameEx <- "abdom"
> 
> ### * abdom
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: abdom
> ### Title: Abdominal Circumference Data
> ### Aliases: abdom
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(abdom)
> attach(abdom)
> plot(x,y)
> detach(abdom)
> 
> 
> 
> cleanEx(); ..nameEx <- "aep"
> 
> ### * aep
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: aep
> ### Title: The Hospital Stay Data
> ### Aliases: aep
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(aep)
> attach(aep)
> pro<-noinap/los
> plot(ward,pro)
> rm(pro)
> detach(aep)
> 
> 
> 
> cleanEx(); ..nameEx <- "aids"
> 
> ### * aids
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: aids
> ### Title: Aids Cases in England and Wales
> ### Aliases: aids
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(aids)
> attach(aids)
> plot(x,y,pch=21,bg=c("red","green3","blue","yellow")[unclass(qrt)])
> detach(aids)
> 
> 
> 
> cleanEx(); ..nameEx <- "bd"
> 
> ### * bd
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: db
> ### Title: Head Circumference of Dutch Boys
> ### Aliases: db
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(db)
> attach(db)
> plot(age,head)
> detach(db)
> 
> 
> 
> cleanEx(); ..nameEx <- "bfp"
> 
> ### * bfp
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: bfp
> ### Title: Functions to fit fractional polynomials in GAMLSS
> ### Aliases: bfp fp pp
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> #fits polynomials with power 1 and .5 
> mod1<-gamlss(y~bfp(x,c(1,0.5)),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4948.285 
GAMLSS-RS iteration 2: Global Deviance = 4948.285 
> # fit the best of one fractional polynomial
> m1<-gamlss(y~fp(x,1),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4967.042 
GAMLSS-RS iteration 2: Global Deviance = 4967.042 
> # fit the best of two fractional polynomials
> m2<-gamlss(y~fp(x,2),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4941.099 
GAMLSS-RS iteration 2: Global Deviance = 4941.099 
> # fit the best of three fractional polynomials
> m3<-gamlss(y~fp(x,3),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4933.937 
GAMLSS-RS iteration 2: Global Deviance = 4933.937 
> 
> 
> 
> cleanEx(); ..nameEx <- "centiles"
> 
> ### * centiles
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: centiles
> ### Title: Plots the centile curves for a GAMLSS object
> ### Aliases: centiles
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> centiles(h,xvar=abdom$x)
% of cases below  0.4 centile is  0.6557377 
% of cases below  2 centile is  2.459016 
% of cases below  10 centile is  8.688525 
% of cases below  25 centile is  26.06557 
% of cases below  50 centile is  50.32787 
% of cases below  75 centile is  74.91803 
% of cases below  90 centile is  88.68852 
% of cases below  98 centile is  98.03279 
% of cases below  99.6 centile is  99.67213 
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "centiles.com"
> 
> ### * centiles.com
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: centiles.com
> ### Title: Comparing centiles from different GAMLSS models
> ### Aliases: centiles.com
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h1<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> h2<-gamlss(y~lo(x,span=0.4), sigma.formula=~lo(x,span=0.4), family=BCT, data=abdom ) 
GAMLSS-RS iteration 1: Global Deviance = 4773.475 
GAMLSS-RS iteration 2: Global Deviance = 4773.18 
GAMLSS-RS iteration 3: Global Deviance = 4773.149 
GAMLSS-RS iteration 4: Global Deviance = 4773.146 
GAMLSS-RS iteration 5: Global Deviance = 4773.146 
> centiles.com(h1,h2,xvar=abdom$x)
********  Model 1 ******** 
% of cases below  0.4 centile is  0.6557377 
% of cases below  10 centile is  8.688525 
% of cases below  50 centile is  50.32787 
% of cases below  90 centile is  88.68852 
% of cases below  99.6 centile is  99.67213 
********  Model 2 ******** 
% of cases below  0.4 centile is  0.3278689 
% of cases below  10 centile is  8.688525 
% of cases below  50 centile is  50.81967 
% of cases below  90 centile is  89.34426 
% of cases below  99.6 centile is  99.83607 
> rm(h1,h2)
> 
> 
> 
> cleanEx(); ..nameEx <- "centiles.pred"
> 
> ### * centiles.pred
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: centiles.pred
> ### Title: Creating predictive centiles values
> ### Aliases: centiles.pred
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> a<-gamlss(y~cs(x),sigma.fo=~cs(x), data=abdom, family=BCT)
GAMLSS-RS iteration 1: Global Deviance = 4772.531 
GAMLSS-RS iteration 2: Global Deviance = 4772.445 
GAMLSS-RS iteration 3: Global Deviance = 4772.429 
GAMLSS-RS iteration 4: Global Deviance = 4772.423 
GAMLSS-RS iteration 5: Global Deviance = 4772.422 
GAMLSS-RS iteration 6: Global Deviance = 4772.422 
> #plot the centiles
> centiles(a,xvar=abdom$x)
% of cases below  0.4 centile is  0.3278689 
% of cases below  2 centile is  2.459016 
% of cases below  10 centile is  8.688525 
% of cases below  25 centile is  26.72131 
% of cases below  50 centile is  50.32787 
% of cases below  75 centile is  74.09836 
% of cases below  90 centile is  90 
% of cases below  98 centile is  98.19672 
% of cases below  99.6 centile is  99.67213 
> # calculate the centiles at new x values 
> newx<-seq(12,40,2)
> mat <- centiles.pred(a, xname="x", xvalues=newx )
> mat
    x      C0.4        C2       C10       C25       C50       C75       C90
1  12  44.32712  47.53252  51.34103  54.21895  57.45386  60.90707  64.39590
2  14  64.82615  69.10822  74.16073  77.95493  82.19681  86.69959  91.22369
3  16  85.80227  90.94380  96.96902 101.46596 106.46692 111.74631 117.02199
4  18 106.90560 112.65921 119.35570 124.32305 129.81801 135.58715 141.32117
5  20 127.75527 133.93173 141.07618 146.34660 152.14918 158.21134 164.20746
6  22 148.08126 154.60111 162.10601 167.61818 173.66415 179.95607 186.15572
7  24 168.12073 175.00558 182.90302 188.68540 195.01077 201.57521 208.02580
8  26 187.58931 194.98697 203.45802 209.65083 216.41621 223.42770 230.30833
9  28 206.01103 214.15915 223.49079 230.31354 237.76788 245.49419 253.07708
10 30 223.17621 232.26930 242.69729 250.33089 258.67980 267.34264 275.85365
11 32 239.35480 249.50728 261.17230 269.72595 279.09479 288.83057 298.40989
12 34 254.84078 266.03010 278.90819 288.36566 298.73791 309.53091 320.16449
13 36 269.69902 281.85532 295.86485 306.16542 317.47378 329.25326 340.87072
14 38 283.90699 296.99350 312.09247 323.20555 335.41677 348.14842 360.71631
15 40 298.22002 312.08705 328.09392 339.88005 352.83540 366.34781 379.69122
         C98     C99.6
1   69.76056  75.10617
2   98.13402 104.96782
3  125.02793 132.88680
4  149.96646 158.39075
5  173.19566 181.89668
6  195.40663 204.31574
7  217.61987 226.82550
8  240.52564 250.31152
9  264.33860 275.12609
10 288.50947 300.64986
11 312.67947 326.39532
12 336.02950 351.30612
13 358.22508 374.95917
14 379.51072 397.65558
15 399.65391 418.93603
> # plot the centiles  
>  mat <- centiles.pred(a, xname="x",xvalues=newx, plot=TRUE )
> # calculate standard-centiles for new values using the fitted model
> newx <- seq(12,40,2)
>  mat <- centiles.pred(a, xname="x",xvalues=newx, type="standard-centiles" )
>  mat
    x        -4        -3        -2        -1         0         1         2
1  12  35.65419  42.31050  47.80738  52.68015  57.45386  62.71870  69.34049
2  14  53.08735  62.11735  69.47411  75.92867  82.19681  89.05186  97.59486
3  16  71.52359  82.53218  91.38161  99.06728 106.46692 114.49285 124.40548
4  18  90.71768 103.22636 113.14743 121.67662 129.81801 138.57609 149.29664
5  20 110.17159 123.78631 134.45419 143.54172 152.14918 161.34045 172.50145
6  22 129.34527 143.87541 155.15124 164.68714 173.66415 183.19428 194.69389
7  24 148.20244 163.66711 175.58548 185.61255 195.01077 204.94663 216.88199
8  26 166.11621 182.79744 195.60951 206.36085 216.41621 227.02501 239.74051
9  28 182.36560 200.73358 214.84489 226.68882 237.76788 249.45856 263.47317
10 30 196.85739 217.29304 233.03510 246.27443 258.67980 271.79114 287.53623
11 32 210.07697 232.79610 250.36312 265.17908 279.09479 293.83572 311.58108
12 34 222.67767 247.62194 266.97416 283.33686 298.73791 315.08520 334.80726
13 36 234.84493 261.86452 282.88167 300.68706 317.47378 335.32000 356.88721
14 38 246.46897 275.48070 298.09903 317.29387 335.41677 354.71009 378.06099
15 40 258.58391 289.29439 313.25880 333.60983 352.83540 373.31380 398.11371
           3         4
1   78.95471  95.37367
2  109.85771 130.46482
3  138.47662 161.75422
4  164.34732 188.86203
5  188.01639 212.94037
6  210.55575 235.76220
7  233.25412 259.07197
8  257.13538 284.46243
9  282.64924 312.78327
10 309.12623 343.15469
11 335.98709 374.61496
12 362.00473 405.21168
13 386.69168 434.17889
14 410.38974 462.03101
15 432.47361 487.41529
> #to plot the centiles  
> mat <- centiles.pred(a, xname="x",xvalues=newx, type="s", plot = TRUE )
> # get z-scores for new values and plot them
> newx <- c(20,21.2,23,20.9,24.2,24.1,25)
> newy <- c(130,121,123,125,140,145,150)
> znewx <- centiles.pred(a, xname="x",xvalues=newx,yval=newy, type="z-scores" )
> znewx
[1] -2.442657 -4.042038 -4.802512 -3.594058 -4.511235 -4.221110 -4.369420
> for(i in 1:7) points(newx[i],newy[i],col="blue")
> # now with transformed data
> aa<-gamlss(y~cs(x^0.5),sigma.fo=~cs(x^0.5), data=abdom, family=BCT)  
GAMLSS-RS iteration 1: Global Deviance = 4771.562 
GAMLSS-RS iteration 2: Global Deviance = 4770.55 
GAMLSS-RS iteration 3: Global Deviance = 4770.524 
GAMLSS-RS iteration 4: Global Deviance = 4770.519 
GAMLSS-RS iteration 5: Global Deviance = 4770.519 
> centiles(aa,xvar=abdom$x)
% of cases below  0.4 centile is  0.3278689 
% of cases below  2 centile is  2.459016 
% of cases below  10 centile is  8.52459 
% of cases below  25 centile is  26.22951 
% of cases below  50 centile is  50.16393 
% of cases below  75 centile is  74.09836 
% of cases below  90 centile is  90.32787 
% of cases below  98 centile is  98.03279 
% of cases below  99.6 centile is  99.67213 
> mat <- centiles.pred(aa, xname="x",xvalues=c(30) )
> xx<-rep(mat[,1],9)
> yy<-mat[,2:10]
> points(xx,yy,col="red")
> # now with transformed previously data
> nx<-abdom$x^0.5
> aa<-gamlss(y~cs(nx),sigma.fo=~cs(nx), data=abdom, family=BCT)
GAMLSS-RS iteration 1: Global Deviance = 4771.562 
GAMLSS-RS iteration 2: Global Deviance = 4770.55 
GAMLSS-RS iteration 3: Global Deviance = 4770.524 
GAMLSS-RS iteration 4: Global Deviance = 4770.519 
GAMLSS-RS iteration 5: Global Deviance = 4770.519 
> centiles(aa, xvar=abdom$x)
% of cases below  0.4 centile is  0.3278689 
% of cases below  2 centile is  2.459016 
% of cases below  10 centile is  8.52459 
% of cases below  25 centile is  26.22951 
% of cases below  50 centile is  50.16393 
% of cases below  75 centile is  74.09836 
% of cases below  90 centile is  90.32787 
% of cases below  98 centile is  98.03279 
% of cases below  99.6 centile is  99.67213 
> newd<-data.frame( abdom, nx=abdom$x^0.5)
> mat <-  centiles.pred(aa, xname="nx", xvalues=c(30), power=0.5, data=newd)
> xxx<-rep(mat[,1],9)
> yyy<-mat[,2:10]
> points(xxx,yyy,col="red")
> 
> 
> 
> cleanEx(); ..nameEx <- "centiles.split"
> 
> ### * centiles.split
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: centiles.split
> ### Title: Plots centile curves split by x for a GAMLSS object
> ### Aliases: centiles.split
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> mout <- centiles.split(h,xvar=abdom$x)
> mout
     12.22 to 20.07 20.07 to 27.07 27.07 to 34.5 34.5 to 42.5
0.4        0.000000       0.000000     0.6410256     2.040816
2          2.597403       1.307190     2.5641026     3.401361
10         9.090909       6.535948     9.6153846     9.523810
25        24.675325      30.065359    23.7179487    25.850340
50        48.051948      53.594771    49.3589744    50.340136
75        73.376623      76.470588    74.3589744    75.510204
90        84.415584      92.156863    88.4615385    89.795918
98        97.402597      99.346405    97.4358974    97.959184
99.6     100.000000      99.346405    99.3589744   100.000000
> rm(h,mout)
> 
> 
> 
> cleanEx(); ..nameEx <- "checklink"
> 
> ### * checklink
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: checklink
> ### Title: Set the Right Link Function for Specified Parameter and
> ###   Distribution
> ### Aliases: checklink
> ### Keywords: regression
> 
> ### ** Examples
> 
> 
> 
> cleanEx(); ..nameEx <- "coef.gamlss"
> 
> ### * coef.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: coef.gamlss
> ### Title: Extract Model Coefficients in a GAMLSS fitted model
> ### Aliases: coef.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> coef(h)
 (Intercept)  poly(x, 3)1  poly(x, 3)2  poly(x, 3)3         qrt2         qrt3 
 4.814472625  8.256092237 -3.353309934  0.938815624 -0.156380604  0.009986146 
        qrt4 
-0.115704890 
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "cs"
> 
> ### * cs
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: cs
> ### Title: Specify a Smoothing Spline Fit in a GAMLSS Formula
> ### Aliases: cs vc
> ### Keywords: regression
> 
> ### ** Examples
> 
> # cubic splines example
> data(aids)
> attach(aids)
> # fitting a smoothing cubic spline with 7 degrees of freedom
> # plus the a quarterly  effect  
> aids1<-gamlss(y~cs(x,df=7)+qrt,data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 374.1748 
GAMLSS-RS iteration 2: Global Deviance = 374.1748 
> plot(x,y)
> lines(x,fitted(aids1))
> rm(aids1)
> detach(aids)
> # varying-coefficient example
> data(rent)
> attach(rent)
> # adjusting the variables
> Flbar<-Fl-mean(Fl)
> Abar<-A-mean(A)
> # additive model
>  m1<-gamlss(R~cs(Flbar, df=3)+cs(Abar))
GAMLSS-RS iteration 1: Global Deviance = 28262.63 
GAMLSS-RS iteration 2: Global Deviance = 28262.63 
> # varying-coefficient model
>  m2<-gamlss(R~cs(Flbar, df=3)+cs(Abar)+vc(r=Abar,x=Flbar))
GAMLSS-RS iteration 1: Global Deviance = 28233.08 
GAMLSS-RS iteration 2: Global Deviance = 28233.08 
> AIC(m1,m2)
         df      AIC
m2 14.00192 28261.08
m1 10.00128 28282.64
> detach(rent)
> 
> 
> 
> cleanEx(); ..nameEx <- "deviance.gamlss"
> 
> ### * deviance.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: deviance.gamlss
> ### Title: Global Deviance of a GAMLSS model
> ### Aliases: deviance.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> deviance(h)
[1] 416.8014
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "fabric"
> 
> ### * fabric
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fabric
> ### Title: The Fabric Data
> ### Aliases: fabric
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(fabric)
> attach(fabric)
> plot(x,y)
> detach(fabric)
> 
> 
> 
> cleanEx(); ..nameEx <- "findhyper"
> 
> ### * findhyper
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: find.hyper
> ### Title: A function to select values of hyperparameters in a GAMLSS model
> ### Aliases: find.hyper
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> attach(abdom)
> # declare the model
> mod1<-quote(gamlss(y~cs(nx,df=p[1]),family=BCT,data=abdom,
+                         control=gamlss.control(trace=FALSE)))
> # we wand also to check for a transformation in x so we use other
> op<-find.hyper(model=mod1, other=quote(nx<-x^p[2]), par=c(3,0.5), 
+             lower=c(1,0.001), steps=c(0.1,0.001))
par 3 0.5 crit= 4814.877 with pen= 2.5 
par 3.1 0.5 crit= 4814.934 with pen= 2.5 
par 2.9 0.5 crit= 4814.834 with pen= 2.5 
par 3 0.501 crit= 4814.88 with pen= 2.5 
par 3 0.499 crit= 4814.875 with pen= 2.5 
par 2.498635 0.001 crit= 4814.813 with pen= 2.5 
par 2.598635 0.001 crit= 4814.592 with pen= 2.5 
par 2.398635 0.001 crit= 4815.131 with pen= 2.5 
par 2.498635 0.002 crit= 4815.198 with pen= 2.5 
par 2.498635 0.001 crit= 4814.813 with pen= 2.5 
par 5.192955 0.001 crit= 4817.336 with pen= 2.5 
par 5.292955 0.001 crit= 4817.501 with pen= 2.5 
par 5.092955 0.001 crit= 4817.172 with pen= 2.5 
par 5.192955 0.002 crit= 4817.373 with pen= 2.5 
par 5.192955 0.001 crit= 4817.338 with pen= 2.5 
par 3.161698 0.001 crit= 4814.394 with pen= 2.5 
par 3.261698 0.001 crit= 4814.466 with pen= 2.5 
par 3.061698 0.001 crit= 4814.343 with pen= 2.5 
par 3.161698 0.002 crit= 4814.679 with pen= 2.5 
par 3.161698 0.001 crit= 4814.393 with pen= 2.5 
par 3.038148 0.001 crit= 4814.335 with pen= 2.5 
par 3.138148 0.001 crit= 4814.38 with pen= 2.5 
par 2.938148 0.001 crit= 4814.319 with pen= 2.5 
par 3.038148 0.002 crit= 4814.638 with pen= 2.5 
par 3.038148 0.001 crit= 4814.334 with pen= 2.5 
par 2.919385 0.001 crit= 4814.32 with pen= 2.5 
par 3.019385 0.001 crit= 4814.33 with pen= 2.5 
par 2.819385 0.001 crit= 4814.348 with pen= 2.5 
par 2.919385 0.002 crit= 4814.64 with pen= 2.5 
par 2.919385 0.001 crit= 4814.319 with pen= 2.5 
par 2.946980 0.001 crit= 4814.32 with pen= 2.5 
par 3.046980 0.001 crit= 4814.338 with pen= 2.5 
par 2.846980 0.001 crit= 4814.336 with pen= 2.5 
par 2.946980 0.002 crit= 4814.635 with pen= 2.5 
par 2.946980 0.001 crit= 4814.318 with pen= 2.5 
par 2.944492 0.001 crit= 4814.319 with pen= 2.5 
par 3.044492 0.001 crit= 4814.337 with pen= 2.5 
par 2.844492 0.001 crit= 4814.337 with pen= 2.5 
par 2.944492 0.002 crit= 4814.635 with pen= 2.5 
par 2.944492 0.001 crit= 4814.318 with pen= 2.5 
par 2.94443 0.001 crit= 4814.320 with pen= 2.5 
par 3.04443 0.001 crit= 4814.337 with pen= 2.5 
par 2.84443 0.001 crit= 4814.337 with pen= 2.5 
par 2.94443 0.002 crit= 4814.635 with pen= 2.5 
par 2.94443 0.001 crit= 4814.318 with pen= 2.5 
par 2.944492 0.001 crit= 4814.320 with pen= 2.5 
par 3.044492 0.001 crit= 4814.337 with pen= 2.5 
par 2.844492 0.001 crit= 4814.337 with pen= 2.5 
par 2.944492 0.002 crit= 4814.635 with pen= 2.5 
par 2.944492 0.001 crit= 4814.318 with pen= 2.5 
par 2.944492 0.001 crit= 4814.320 with pen= 2.5 
par 3.044492 0.001 crit= 4814.337 with pen= 2.5 
par 2.844492 0.001 crit= 4814.337 with pen= 2.5 
par 2.944492 0.002 crit= 4814.635 with pen= 2.5 
par 2.944492 0.001 crit= 4814.318 with pen= 2.5 
par 2.944492 0.001 crit= 4814.320 with pen= 2.5 
par 3.044492 0.001 crit= 4814.337 with pen= 2.5 
par 2.844492 0.001 crit= 4814.337 with pen= 2.5 
par 2.944492 0.002 crit= 4814.635 with pen= 2.5 
par 2.944492 0.001 crit= 4814.318 with pen= 2.5 
par 2.944492 0.001 crit= 4814.320 with pen= 2.5 
par 3.044492 0.001 crit= 4814.337 with pen= 2.5 
par 2.844492 0.001 crit= 4814.337 with pen= 2.5 
par 2.944492 0.002 crit= 4814.635 with pen= 2.5 
par 2.944492 0.001 crit= 4814.318 with pen= 2.5 
> # the optimum parameters found are 
> # par=(p[1],p[2]) = (2.944836 0.001000) = (df for mu, lambda)
> # so it needs df = 3 on top of the constant and linear 
> # and an log transformation for x
> op
$par
[1] 2.944492 0.001000

$value
[1] 4814.320

$counts
function gradient 
      13       13 

$convergence
[1] 0

$message
[1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

> rm(op)
> 
> 
> 
> cleanEx(); ..nameEx <- "fitted.gamlss"
> 
> ### * fitted.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fitted.gamlss
> ### Title: Extract Fitted Values For A GAMLSS Model
> ### Aliases: fitted.gamlss fv
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> fitted(h)
 [1]   3.896026   4.614935   7.422221   8.768654  12.979298  14.411163
 [7]  21.768152  24.200524  33.775141  35.428204  50.655403  53.411295
[13]  70.836914  70.748250  96.504249  97.265284 123.548781 118.412788
[19] 155.304275 150.799148 184.898568 171.394836 217.839470 205.379020
[25] 244.987968 221.367038 274.793082 253.529967 296.532272 263.235310
[31] 321.655483 292.697106 338.310168 297.365910 360.488934 326.079382
[37] 375.382385 329.271344 399.124699 361.696170 417.973122 368.750221
[43] 450.444031 412.173451 481.879358
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "fitted.plot"
> 
> ### * fitted.plot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: fitted.plot
> ### Title: Plots The Fitted Values of a GAMLSS Model
> ### Aliases: fitted.plot
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h1<-gamlss(y~cs(x,df=3), sigma.formula=~x, family=BCT, data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4789.114 
GAMLSS-RS iteration 2: Global Deviance = 4788.652 
GAMLSS-RS iteration 3: Global Deviance = 4788.64 
GAMLSS-RS iteration 4: Global Deviance = 4788.639 
GAMLSS-RS iteration 5: Global Deviance = 4788.639 
> h2<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,df=3), family=BCT, data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4772.531 
GAMLSS-RS iteration 2: Global Deviance = 4772.445 
GAMLSS-RS iteration 3: Global Deviance = 4772.429 
GAMLSS-RS iteration 4: Global Deviance = 4772.423 
GAMLSS-RS iteration 5: Global Deviance = 4772.422 
GAMLSS-RS iteration 6: Global Deviance = 4772.422 
> fitted.plot(h1,h2,x=abdom$x)
> rm(h1,h2)
> 
> 
> 
> cleanEx(); ..nameEx <- "formula.gamlss"
> 
> ### * formula.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: formula.gamlss
> ### Title: Extract the Model Formula in a GAMLSS fitted model
> ### Aliases: formula.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> formula(h,"mu")
y ~ poly(x, 3) + qrt
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "gamlss"
> 
> ### * gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gamlss
> ### Title: Generalized Additive Models for Location Scale and Shape
> ### Aliases: gamlss is.gamlss gamlssNews
> ### Keywords: regression
> 
> ### ** Examples
>  
> data(abdom)
> mod<-gamlss(y~cs(x,df=3),sigma.fo=~cs(x,df=1),family=BCT, data=abdom, method=mixed(1,20))
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-CG iteration 1: Global Deviance = 4775.962 
GAMLSS-CG iteration 2: Global Deviance = 4775.876 
GAMLSS-CG iteration 3: Global Deviance = 4775.877 
> plot(mod)
*******************************************************************
	      Summary of the Quantile Residuals
                           mean   =  0.001107079 
                       variance   =  1.001849 
               coef. of skewness  =  -0.009695628 
               coef. of kurtosis  =  2.992983 
Filliben correlation coefficient  =  0.9992747 
*******************************************************************
> rm(mod)
> 
> 
> 
> cleanEx(); ..nameEx <- "gamlss.control"
> 
> ### * gamlss.control
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gamlss.control
> ### Title: Auxiliary for Controlling GAMLSS Fitting
> ### Aliases: gamlss.control
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> con<-gamlss.control(mu.step=0.1)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids, control=con) # 
GAMLSS-RS iteration 1: Global Deviance = 416.835 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
GAMLSS-RS iteration 3: Global Deviance = 416.8014 
> rm(h,con)
> 
> 
> 
> cleanEx(); ..nameEx <- "gamlss.family"
> 
> ### * gamlss.family
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gamlss.family
> ### Title: Family Objects for fitting a GAMLSS model
> ### Aliases: gamlss.family as.gamlss.family print.gamlss.family
> ###   gamlss.family.default as.family
> ### Keywords: regression
> 
> ### ** Examples
> 
>  normal<-NO(mu.link="log", sigma.link="log")
>  normal

GAMLSS Family: NO Normal 
Link function for mu   : log 
Link function for sigma: log 
>  rm(normal)
> 
> 
> 
> cleanEx(); ..nameEx <- "gamlss.scope"
> 
> ### * gamlss.scope
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: gamlss.scope
> ### Title: Generate a Scope Argument for Stepwise GAMLSS
> ### Aliases: gamlss.scope
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(usair)
> gs1<-gamlss.scope(model.frame(y~x1+x2+x3+x4+x5+x6, data=usair))
> gs2<-gamlss.scope(model.frame(usair))
> gs1
$x1
~1 + x1 + cs(x1)

$x2
~1 + x2 + cs(x2)

$x3
~1 + x3 + cs(x3)

$x4
~1 + x4 + cs(x4)

$x5
~1 + x5 + cs(x5)

$x6
~1 + x6 + cs(x6)

> gs2
$x1
~1 + x1 + cs(x1)

$x2
~1 + x2 + cs(x2)

$x3
~1 + x3 + cs(x3)

$x4
~1 + x4 + cs(x4)

$x5
~1 + x5 + cs(x5)

$x6
~1 + x6 + cs(x6)

> gs3<-gamlss.scope(model.frame(usair), smooth="fp", arg="3")
> gs3 
$x1
~1 + x1 + fp(x1, 3)

$x2
~1 + x2 + fp(x2, 3)

$x3
~1 + x3 + fp(x3, 3)

$x4
~1 + x4 + fp(x4, 3)

$x5
~1 + x5 + fp(x5, 3)

$x6
~1 + x6 + fp(x6, 3)

> 
> 
> 
> cleanEx(); ..nameEx <- "glim.control"
> 
> ### * glim.control
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: glim.control
> ### Title: Auxiliary for Controlling the inner algorithm in a GAMLSS
> ###   Fitting
> ### Aliases: glim.control
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> con<-glim.control(trace=TRUE)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids, i.control=con) # 
GLIM iteration 1: Global Deviance = 785.86 
GLIM iteration 2: Global Deviance = 473.255 
GLIM iteration 3: Global Deviance = 419.3565 
GLIM iteration 4: Global Deviance = 416.8088 
GLIM iteration 5: Global Deviance = 416.8014 
GLIM iteration 6: Global Deviance = 416.8014 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GLIM iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> rm(h,con)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "histDist"
> 
> ### * histDist
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: histDist
> ### Title: This function plots the histogram and a fitted (GAMLSS family)
> ###   distribution to a variable
> ### Aliases: histDist
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> attach(abdom)
> histDist(y,family="NO")
GAMLSS-RS iteration 1: Global Deviance = 7201.417 
GAMLSS-RS iteration 2: Global Deviance = 7201.417 
> # use the ymax argument of the of the truehist()
> histDist(y,family="NO",ymax=0.005)
GAMLSS-RS iteration 1: Global Deviance = 7201.417 
GAMLSS-RS iteration 2: Global Deviance = 7201.417 
> # bad fit use PE
> histDist(y,family="PE",ymax=0.005)
GAMLSS-RS iteration 1: Global Deviance = 7093.59 
GAMLSS-RS iteration 2: Global Deviance = 7093.034 
GAMLSS-RS iteration 3: Global Deviance = 7093.013 
GAMLSS-RS iteration 4: Global Deviance = 7093.007 
GAMLSS-RS iteration 5: Global Deviance = 7093.005 
GAMLSS-RS iteration 6: Global Deviance = 7093.004 
> detach(abdom)
> 
> 
> 
> cleanEx(); ..nameEx <- "hodges"
> 
> ### * hodges
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: hodges
> ### Title: Hodges data
> ### Aliases: hodges hodges1
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(hodges)
> attach(hodges)
> plot(prind~state, cex=1, cex.lab=1.5, cex.axis=1, cex.main=1.2)
> str(hodges)
`data.frame':	341 obs. of  6 variables:
 $ state: Factor w/ 45 levels "AL","AZ","CA",..: 1 1 1 1 1 1 1 2 2 2 ...
 $ plan : Factor w/ 325 levels "11","12","15",..: 39 77 216 218 221 223 237 4 37 38 ...
 $ prind: num  236 168 170 181 164 ...
 $ prfam: num  470 429 411 424 417 ...
 $ enind: num   93 824 483 174 285 ...
 $ enfam: num   321 1508  789  338  438 ...
> data(hodges1)
> str(hodges1)
`data.frame':	45 obs. of  4 variables:
 $ State : Factor w/ 45 levels "AL","AZ","CA",..: 1 2 3 4 5 6 7 8 9 10 ...
 $ expe  : num  4455 5212 6169 5449 6696 ...
 $ pop   : num   4041  3665 29760  3294  3287 ...
 $ region: Factor w/ 7 levels "MA","MT","NC",..: 7 2 5 2 4 6 6 6 6 5 ...
> 
> 
> 
> cleanEx(); ..nameEx <- "lo"
> 
> ### * lo
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: lo
> ### Title: Specify a loess fit in a GAMLSS formula
> ### Aliases: lo
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> attach(aids)
> # fitting a loess curve with span=0.4 plus the a quarterly  effect 
> aids1<-gamlss(y~lo(x,span=0.4)+qrt,data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 399.7146 
GAMLSS-RS iteration 2: Global Deviance = 399.747 
GAMLSS-RS iteration 3: Global Deviance = 399.7476 
> plot(x,y)
> lines(x,fitted(aids1))
> rm(aids1)
> detach(aids)
> 
> 
> 
> cleanEx(); ..nameEx <- "lpred"
> 
> ### * lpred
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: lpred
> ### Title: Extract Linear Predictor Values and Standard Errors For A GAMLSS
> ###   Model
> ### Aliases: lpred lp
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> mod<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> mod.t <- lpred(mod, type = "terms", terms= "qrt")
> mod.t
           qrt
1   0.06406873
2  -0.09231187
3   0.07405488
4  -0.05163616
5   0.06406873
6  -0.09231187
7   0.07405488
8  -0.05163616
9   0.06406873
10 -0.09231187
11  0.07405488
12 -0.05163616
13  0.06406873
14 -0.09231187
15  0.07405488
16 -0.05163616
17  0.06406873
18 -0.09231187
19  0.07405488
20 -0.05163616
21  0.06406873
22 -0.09231187
23  0.07405488
24 -0.05163616
25  0.06406873
26 -0.09231187
27  0.07405488
28 -0.05163616
29  0.06406873
30 -0.09231187
31  0.07405488
32 -0.05163616
33  0.06406873
34 -0.09231187
35  0.07405488
36 -0.05163616
37  0.06406873
38 -0.09231187
39  0.07405488
40 -0.05163616
41  0.06406873
42 -0.09231187
43  0.07405488
44 -0.05163616
45  0.06406873
attr(,"constant")
[1] 4.750404
> mod.lp <- lp(mod)
> mod.lp 
       1        2        3        4        5        6        7        8 
1.359957 1.529298 2.004478 2.171183 2.563356 2.668003 3.080448 3.186374 
       9       10       11       12       13       14       15       16 
3.519725 3.567508 3.925046 3.978022 4.260380 4.259128 4.569587 4.577442 
      17       18       19       20       21       22       23       24 
4.816636 4.774177 5.045386 5.015949 5.219807 5.143970 5.383758 5.324857 
      25       26       27       28       29       30       31       32 
5.501209 5.399822 5.616018 5.535482 5.692156 5.573048 5.773481 5.679138 
      33       34       35       36       37       38       39       40 
5.823963 5.694963 5.887461 5.787141 5.927945 5.796882 5.989274 5.890805 
      41       42       43       44       45 
6.035417 5.910120 6.110234 6.021444 6.177694 
> rm(mod, mod.t,mod.lp)
> 
> 
> 
> cleanEx(); ..nameEx <- "make.link.gamlss"
> 
> ### * make.link.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: make.link.gamlss
> ### Title: Create a Link for GAMLSS families
> ### Aliases: make.link.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> str(make.link.gamlss("logitshifted"), par=c(-1,1))
List of 4
 $ linkfun :function (mu, shift = par)  
 $ linkinv :function (eta, shift = par)  
 $ mu.eta  :function (eta, shift = par)  
 $ valideta:function (eta)  
> l2<-make.link.gamlss("logitshifted", par=c(-1,1))
> l2$linkfun(0) # should be zero
[1] 0
> l2$linkfun(1) # should be Inf
[1] Inf
> l2$linkinv(-5:5) 
 [1] -0.9866143 -0.9640276 -0.9051483 -0.7615942 -0.4621172  0.0000000
 [7]  0.4621172  0.7615942  0.9051483  0.9640276  0.9866143
> 
> 
> 
> cleanEx(); ..nameEx <- "model.frame.gamlss"
> 
> ### * model.frame.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: model.frame.gamlss
> ### Title: Extract a model.frame, a model matrix or terms from a GAMLSS
> ###   object for a given distributional parameter
> ### Aliases: model.frame.gamlss terms.gamlss model.matrix.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> mod<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> model.frame(mod)
     y  poly(x, 3).1  poly(x, 3).2  poly(x, 3).3 qrt
1    2 -2.525235e-01  3.118254e-01 -3.451275e-01   1
2    6 -2.410452e-01  2.693037e-01 -2.510018e-01   2
3   10 -2.295668e-01  2.287598e-01 -1.678210e-01   3
4    8 -2.180885e-01  1.901937e-01 -9.506380e-02   4
5   12 -2.066101e-01  1.536053e-01 -3.220912e-02   1
6    9 -1.951318e-01  1.189947e-01  2.126427e-02   2
7   28 -1.836535e-01  8.636179e-02  6.587755e-02   3
8   28 -1.721751e-01  5.570665e-02  1.021519e-01   4
9   35 -1.606968e-01  2.702926e-02  1.306085e-01   1
10  32 -1.492184e-01  3.296251e-04  1.517685e-01   2
11  46 -1.377401e-01 -2.439226e-02  1.661532e-01   3
12  47 -1.262617e-01 -4.713640e-02  1.742836e-01   4
13  50 -1.147834e-01 -6.790278e-02  1.766811e-01   1
14  61 -1.033051e-01 -8.669141e-02  1.738667e-01   2
15  99 -9.182673e-02 -1.035023e-01  1.663616e-01   3
16  95 -8.034839e-02 -1.183354e-01  1.546871e-01   4
17 150 -6.887004e-02 -1.311908e-01  1.393644e-01   1
18 143 -5.739170e-02 -1.420684e-01  1.209145e-01   2
19 197 -4.591336e-02 -1.509683e-01  9.985868e-02   3
20 159 -3.443502e-02 -1.578904e-01  7.671815e-02   4
21 204 -2.295668e-02 -1.628348e-01  5.201407e-02   1
22 168 -1.147834e-02 -1.658014e-01  2.626763e-02   2
23 196 -2.075173e-18 -1.667903e-01 -4.048628e-17   3
24 194  1.147834e-02 -1.658014e-01 -2.626763e-02   4
25 210  2.295668e-02 -1.628348e-01 -5.201407e-02   1
26 180  3.443502e-02 -1.578904e-01 -7.671815e-02   2
27 277  4.591336e-02 -1.509683e-01 -9.985868e-02   3
28 181  5.739170e-02 -1.420684e-01 -1.209145e-01   4
29 327  6.887004e-02 -1.311908e-01 -1.393644e-01   1
30 276  8.034839e-02 -1.183354e-01 -1.546871e-01   2
31 365  9.182673e-02 -1.035023e-01 -1.663616e-01   3
32 300  1.033051e-01 -8.669141e-02 -1.738667e-01   4
33 356  1.147834e-01 -6.790278e-02 -1.766811e-01   1
34 304  1.262617e-01 -4.713640e-02 -1.742836e-01   2
35 307  1.377401e-01 -2.439226e-02 -1.661532e-01   3
36 386  1.492184e-01  3.296251e-04 -1.517685e-01   4
37 331  1.606968e-01  2.702926e-02 -1.306085e-01   1
38 358  1.721751e-01  5.570665e-02 -1.021519e-01   2
39 415  1.836535e-01  8.636179e-02 -6.587755e-02   3
40 374  1.951318e-01  1.189947e-01 -2.126427e-02   4
41 412  2.066101e-01  1.536053e-01  3.220912e-02   1
42 358  2.180885e-01  1.901937e-01  9.506380e-02   2
43 416  2.295668e-01  2.287598e-01  1.678210e-01   3
44 414  2.410452e-01  2.693037e-01  2.510018e-01   4
45 496  2.525235e-01  3.118254e-01  3.451275e-01   1
> model.matrix(mod)
   (Intercept)   poly(x, 3)1   poly(x, 3)2   poly(x, 3)3 qrt2 qrt3 qrt4
1            1 -2.525235e-01  0.3118253845 -3.451275e-01    0    0    0
2            1 -2.410452e-01  0.2693037411 -2.510018e-01    1    0    0
3            1 -2.295668e-01  0.2287598487 -1.678210e-01    0    1    0
4            1 -2.180885e-01  0.1901937070 -9.506380e-02    0    0    1
5            1 -2.066101e-01  0.1536053162 -3.220912e-02    0    0    0
6            1 -1.951318e-01  0.1189946763  2.126427e-02    1    0    0
7            1 -1.836535e-01  0.0863617872  6.587755e-02    0    1    0
8            1 -1.721751e-01  0.0557066490  1.021519e-01    0    0    1
9            1 -1.606968e-01  0.0270292617  1.306085e-01    0    0    0
10           1 -1.492184e-01  0.0003296251  1.517685e-01    1    0    0
11           1 -1.377401e-01 -0.0243922605  1.661532e-01    0    1    0
12           1 -1.262617e-01 -0.0471363953  1.742836e-01    0    0    1
13           1 -1.147834e-01 -0.0679027793  1.766811e-01    0    0    0
14           1 -1.033051e-01 -0.0866914124  1.738667e-01    1    0    0
15           1 -9.182673e-02 -0.1035022946  1.663616e-01    0    1    0
16           1 -8.034839e-02 -0.1183354260  1.546871e-01    0    0    1
17           1 -6.887004e-02 -0.1311908066  1.393644e-01    0    0    0
18           1 -5.739170e-02 -0.1420684363  1.209145e-01    1    0    0
19           1 -4.591336e-02 -0.1509683151  9.985868e-02    0    1    0
20           1 -3.443502e-02 -0.1578904431  7.671815e-02    0    0    1
21           1 -2.295668e-02 -0.1628348202  5.201407e-02    0    0    0
22           1 -1.147834e-02 -0.1658014465  2.626763e-02    1    0    0
23           1 -2.075173e-18 -0.1667903219 -4.048628e-17    0    1    0
24           1  1.147834e-02 -0.1658014465 -2.626763e-02    0    0    1
25           1  2.295668e-02 -0.1628348202 -5.201407e-02    0    0    0
26           1  3.443502e-02 -0.1578904431 -7.671815e-02    1    0    0
27           1  4.591336e-02 -0.1509683151 -9.985868e-02    0    1    0
28           1  5.739170e-02 -0.1420684363 -1.209145e-01    0    0    1
29           1  6.887004e-02 -0.1311908066 -1.393644e-01    0    0    0
30           1  8.034839e-02 -0.1183354260 -1.546871e-01    1    0    0
31           1  9.182673e-02 -0.1035022946 -1.663616e-01    0    1    0
32           1  1.033051e-01 -0.0866914124 -1.738667e-01    0    0    1
33           1  1.147834e-01 -0.0679027793 -1.766811e-01    0    0    0
34           1  1.262617e-01 -0.0471363953 -1.742836e-01    1    0    0
35           1  1.377401e-01 -0.0243922605 -1.661532e-01    0    1    0
36           1  1.492184e-01  0.0003296251 -1.517685e-01    0    0    1
37           1  1.606968e-01  0.0270292617 -1.306085e-01    0    0    0
38           1  1.721751e-01  0.0557066490 -1.021519e-01    1    0    0
39           1  1.836535e-01  0.0863617872 -6.587755e-02    0    1    0
40           1  1.951318e-01  0.1189946763 -2.126427e-02    0    0    1
41           1  2.066101e-01  0.1536053162  3.220912e-02    0    0    0
42           1  2.180885e-01  0.1901937070  9.506380e-02    1    0    0
43           1  2.295668e-01  0.2287598487  1.678210e-01    0    1    0
44           1  2.410452e-01  0.2693037411  2.510018e-01    0    0    1
45           1  2.525235e-01  0.3118253845  3.451275e-01    0    0    0
attr(,"assign")
[1] 0 1 1 1 2 2 2
attr(,"contrasts")
attr(,"contrasts")$qrt
[1] "contr.treatment"

> terms(mod, "mu")
y ~ poly(x, 3) + qrt
attr(,"variables")
list(y, poly(x, 3), qrt)
attr(,"factors")
           poly(x, 3) qrt
y                   0   0
poly(x, 3)          1   0
qrt                 0   1
attr(,"term.labels")
[1] "poly(x, 3)" "qrt"       
attr(,"specials")
attr(,"specials")$cs
NULL

attr(,"specials")$lo
NULL

attr(,"specials")$random
NULL

attr(,"specials")$ra
NULL

attr(,"specials")$rc
NULL

attr(,"specials")$fp
NULL

attr(,"specials")$ps
NULL

attr(,"specials")$vc
NULL

attr(,"specials")$pp
NULL

attr(,"order")
[1] 1 1
attr(,"intercept")
[1] 1
attr(,"response")
[1] 1
attr(,".Environment")
<environment: R_GlobalEnv>
attr(,"predvars")
list(y, poly(x, 3, coefs = list(alpha = c(23, 23, 23), norm2 = c(1, 
45, 7590, 1022626, 132532329.6))), qrt)
attr(,"dataClasses")
          y  poly(x, 3)         qrt 
  "numeric" "nmatrix.3"    "factor" 
> rm(mod)
> 
> 
> 
> cleanEx(); ..nameEx <- "par.plot"
> 
> ### * par.plot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: par.plot
> ### Title: A function to plot parallel plot for repeated measurement data
> ### Aliases: par.plot
> ### Keywords: regression
> 
> ### ** Examples
> 
> 
> library(nlme)
> data(Orthodont)
> par.plot(distance~age,data=Orthodont,sub=Subject)
> par.plot(distance~age|Sex,data=Orthodont,sub=Subject)
> par.plot(distance~age|Subject,data=Orthodont,sub=Subject,show.given=FALSE)
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "pdf.plot"
> 
> ### * pdf.plot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: pdf.plot
> ### Title: Plots Probability Distribution Functions for GAMLSS Family
> ### Aliases: pdf.plot
> ### Keywords: regression
> 
> ### ** Examples
> 
> pdf.plot(family=BCT, min=1, max=20, step=.05, mu=10, sigma=0.15, nu=-1, tau=c(4,10,20,40) )
> # now using an gamlss object
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) # fits 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> pdf.plot(obj=h , obs=c(23,67), min=50, max=150, step=.5)
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.gamlss"
> 
> ### * plot.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.gamlss
> ### Title: Plot Residual Diagnostics for an GAMLSS Object
> ### Aliases: plot.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> 
> data(aids)
> a<-gamlss(y~cs(x,df=7)+qrt,family=PO,data=aids)
GAMLSS-RS iteration 1: Global Deviance = 374.1748 
GAMLSS-RS iteration 2: Global Deviance = 374.1748 
> plot(a)
*******************************************************************
	 Summary of the Randomised Quantile Residuals
                           mean   =  -0.02949705 
                       variance   =  1.776076 
               coef. of skewness  =  -0.5494088 
               coef. of kurtosis  =  3.766308 
Filliben correlation coefficient  =  0.9831796 
*******************************************************************
> rm(a)
> 
> 
> 
> cleanEx(); ..nameEx <- "predict.gamlss"
> 
> ### * predict.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: predict.gamlss
> ### Title: Extract Predictor Values and Standard Errors For New Data In A
> ###   GAMLSS Model
> ### Aliases: predict.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> a<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> newaids<-data.frame(x=c(45,46,47), qrt=c(2,3,4))
> ap <- predict(a, newdata=newaids, type = "response")
> ap
[1] 412.1194 508.9537 471.5216
> rm(a, ap)
> data(abdom)
> # transform x 
> aa<-gamlss(y~cs(x^.5),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4936.53 
GAMLSS-RS iteration 2: Global Deviance = 4936.53 
> # predict at old values
> predict(aa)[610]
[1] 371.3929
> # predict at new values 
> predict(aa,newdata=data.frame(x=42.43))
[1] 371.3929
> # now transform x first 
> nx<-abdom$x^.5
> aaa<-gamlss(y~cs(nx),data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4936.53 
GAMLSS-RS iteration 2: Global Deviance = 4936.53 
> # create a new data frame 
> newd<-data.frame( abdom, nx=abdom$x^0.5)
> # predict at old values
> predict(aaa)[610]
[1] 371.3929
> # predict at new values 
> predict(aaa,newdata=data.frame(nx=42.43^.5), data=newd)
[1] 371.3929
> 
> 
> 
> cleanEx(); ..nameEx <- "print.gamlss"
> 
> ### * print.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: print.gamlss
> ### Title: Prints a GAMLSS fitted model
> ### Aliases: print.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> print(h) # or just h

Family:  c("PO", "Poisson") 
Fitting method: RS() 

Call:  gamlss(formula = y ~ poly(x, 3) + qrt, family = PO, data = aids) 

Mu Coefficients:
(Intercept)  poly(x, 3)1  poly(x, 3)2  poly(x, 3)3         qrt2         qrt3  
   4.814473     8.256092    -3.353310     0.938816    -0.156381     0.009986  
       qrt4  
  -0.115705  

Degrees of Freedom for the fit: 7 Residual Deg. of Freedom   38 
Global Deviance:     416.801 
            AIC:     430.801 
            SBC:     443.448 
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "prof.dev"
> 
> ### * prof.dev
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: prof.dev
> ### Title: Plotting the Profile Deviance for one of the Parameters in a
> ###   GAMLSS model
> ### Aliases: prof.dev
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> h<-gamlss(y~cs(x,df=3), sigma.formula=~cs(x,1), family=BCT, data=abdom) 
GAMLSS-RS iteration 1: Global Deviance = 4776.174 
GAMLSS-RS iteration 2: Global Deviance = 4775.895 
GAMLSS-RS iteration 3: Global Deviance = 4775.876 
GAMLSS-RS iteration 4: Global Deviance = 4775.876 
> prof.dev(h,"nu",min=-2.000,max=2,step=0.25,type="l")
******************************************************************* 
nu.start=( -2 ) 
GAMLSS-RS iteration 1: Global Deviance = 4809.268 
GAMLSS-RS iteration 2: Global Deviance = 4784.009 
GAMLSS-RS iteration 3: Global Deviance = 4784.094 
GAMLSS-RS iteration 4: Global Deviance = 4784.105 
GAMLSS-RS iteration 5: Global Deviance = 4784.106 
GAMLSS-RS iteration 6: Global Deviance = 4784.107 
******************************************************************* 
nu.start=( -1.75 ) 
GAMLSS-RS iteration 1: Global Deviance = 4782.149 
GAMLSS-RS iteration 2: Global Deviance = 4782.139 
GAMLSS-RS iteration 3: Global Deviance = 4782.150 
GAMLSS-RS iteration 4: Global Deviance = 4782.152 
GAMLSS-RS iteration 5: Global Deviance = 4782.153 
GAMLSS-RS iteration 6: Global Deviance = 4782.154 
******************************************************************* 
nu.start=( -1.5 ) 
GAMLSS-RS iteration 1: Global Deviance = 4780.429 
GAMLSS-RS iteration 2: Global Deviance = 4780.426 
GAMLSS-RS iteration 3: Global Deviance = 4780.43 
GAMLSS-RS iteration 4: Global Deviance = 4780.432 
GAMLSS-RS iteration 5: Global Deviance = 4780.433 
GAMLSS-RS iteration 6: Global Deviance = 4780.433 
******************************************************************* 
nu.start=( -1.25 ) 
GAMLSS-RS iteration 1: Global Deviance = 4778.969 
GAMLSS-RS iteration 2: Global Deviance = 4778.966 
GAMLSS-RS iteration 3: Global Deviance = 4778.970 
GAMLSS-RS iteration 4: Global Deviance = 4778.971 
GAMLSS-RS iteration 5: Global Deviance = 4778.972 
GAMLSS-RS iteration 6: Global Deviance = 4778.973 
******************************************************************* 
nu.start=( -1 ) 
GAMLSS-RS iteration 1: Global Deviance = 4777.776 
GAMLSS-RS iteration 2: Global Deviance = 4777.774 
GAMLSS-RS iteration 3: Global Deviance = 4777.777 
GAMLSS-RS iteration 4: Global Deviance = 4777.779 
GAMLSS-RS iteration 5: Global Deviance = 4777.78 
******************************************************************* 
nu.start=( -0.75 ) 
GAMLSS-RS iteration 1: Global Deviance = 4776.863 
GAMLSS-RS iteration 2: Global Deviance = 4776.862 
GAMLSS-RS iteration 3: Global Deviance = 4776.865 
GAMLSS-RS iteration 4: Global Deviance = 4776.866 
GAMLSS-RS iteration 5: Global Deviance = 4776.867 
******************************************************************* 
nu.start=( -0.5 ) 
GAMLSS-RS iteration 1: Global Deviance = 4776.241 
GAMLSS-RS iteration 2: Global Deviance = 4776.241 
******************************************************************* 
nu.start=( -0.25 ) 
GAMLSS-RS iteration 1: Global Deviance = 4775.918 
GAMLSS-RS iteration 2: Global Deviance = 4775.919 
******************************************************************* 
nu.start=( 1e-04 ) 
GAMLSS-RS iteration 1: Global Deviance = 4775.904 
GAMLSS-RS iteration 2: Global Deviance = 4775.905 
GAMLSS-RS iteration 3: Global Deviance = 4775.907 
GAMLSS-RS iteration 4: Global Deviance = 4775.908 
******************************************************************* 
nu.start=( 0.25 ) 
GAMLSS-RS iteration 1: Global Deviance = 4776.203 
GAMLSS-RS iteration 2: Global Deviance = 4776.205 
GAMLSS-RS iteration 3: Global Deviance = 4776.205 
******************************************************************* 
nu.start=( 0.5 ) 
GAMLSS-RS iteration 1: Global Deviance = 4776.816 
GAMLSS-RS iteration 2: Global Deviance = 4776.816 
******************************************************************* 
nu.start=( 0.75 ) 
GAMLSS-RS iteration 1: Global Deviance = 4777.741 
GAMLSS-RS iteration 2: Global Deviance = 4777.739 
GAMLSS-RS iteration 3: Global Deviance = 4777.739 
******************************************************************* 
nu.start=( 1 ) 
GAMLSS-RS iteration 1: Global Deviance = 4778.973 
GAMLSS-RS iteration 2: Global Deviance = 4778.969 
GAMLSS-RS iteration 3: Global Deviance = 4778.97 
******************************************************************* 
nu.start=( 1.25 ) 
GAMLSS-RS iteration 1: Global Deviance = 4780.504 
GAMLSS-RS iteration 2: Global Deviance = 4780.499 
GAMLSS-RS iteration 3: Global Deviance = 4780.499 
******************************************************************* 
nu.start=( 1.5 ) 
GAMLSS-RS iteration 1: Global Deviance = 4782.326 
GAMLSS-RS iteration 2: Global Deviance = 4782.324 
GAMLSS-RS iteration 3: Global Deviance = 4782.316 
GAMLSS-RS iteration 4: Global Deviance = 4782.317 
******************************************************************* 
nu.start=( 1.75 ) 
GAMLSS-RS iteration 1: Global Deviance = 4784.425 
GAMLSS-RS iteration 2: Global Deviance = 4784.423 
GAMLSS-RS iteration 3: Global Deviance = 4784.414 
GAMLSS-RS iteration 4: Global Deviance = 4784.415 
GAMLSS-RS iteration 5: Global Deviance = 4784.413 
GAMLSS-RS iteration 6: Global Deviance = 4784.414 
******************************************************************* 
nu.start=( 2 ) 
GAMLSS-RS iteration 1: Global Deviance = 4786.787 
GAMLSS-RS iteration 2: Global Deviance = 4786.785 
GAMLSS-RS iteration 3: Global Deviance = 4786.775 
GAMLSS-RS iteration 4: Global Deviance = 4786.777 
GAMLSS-RS iteration 5: Global Deviance = 4786.775 
GAMLSS-RS iteration 6: Global Deviance = 4786.775 
******************************************************************* 
******************************************************************* 
Best estimate of the fixed parameter is  1e-04 
with a Global Deviance equal to  4775.908  at position  9 
A  95 %  Confidence interval is: ( -1.38295 , 1.127437 ) 
******************************************************************* 
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "prof.term"
> 
> ### * prof.term
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: prof.term
> ### Title: Plotting the Profile: deviance or information criterion for one
> ###   of the terms (or hyper-parameters) in a GAMLSS model
> ### Aliases: prof.term
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> gamlss(y~x+qrt,family=NBI,data=aids)
GAMLSS-RS iteration 1: Global Deviance = 492.7247 
GAMLSS-RS iteration 2: Global Deviance = 492.6375 
GAMLSS-RS iteration 3: Global Deviance = 492.6373 

Family:  c("NBI", "Negative Binomial type I") 
Fitting method: RS() 

Call:  gamlss(formula = y ~ x + qrt, family = NBI, data = aids) 

Mu Coefficients:
(Intercept)            x         qrt2         qrt3         qrt4  
    2.88540      0.08744     -0.12038      0.11176     -0.07554  
Sigma Coefficients:
(Intercept)  
     -1.603  

Degrees of Freedom for the fit: 6 Residual Deg. of Freedom   39 
Global Deviance:     492.637 
            AIC:     504.637 
            SBC:     515.477 
> mod<-quote(gamlss(y ~ offset(this * x) + qrt, data = aids, family = NBI))
> prof.term(mod, min=0.06, max=0.11, step=0.001)
GAMLSS-RS iteration 1: Global Deviance = 508.1867 
GAMLSS-RS iteration 2: Global Deviance = 508.1845 
GAMLSS-RS iteration 3: Global Deviance = 508.1845 
GAMLSS-RS iteration 1: Global Deviance = 507.2312 
GAMLSS-RS iteration 2: Global Deviance = 507.2312 
GAMLSS-RS iteration 1: Global Deviance = 506.293 
GAMLSS-RS iteration 2: Global Deviance = 506.2929 
GAMLSS-RS iteration 1: Global Deviance = 505.3713 
GAMLSS-RS iteration 2: Global Deviance = 505.3712 
GAMLSS-RS iteration 1: Global Deviance = 504.4678 
GAMLSS-RS iteration 2: Global Deviance = 504.4677 
GAMLSS-RS iteration 1: Global Deviance = 503.5842 
GAMLSS-RS iteration 2: Global Deviance = 503.5841 
GAMLSS-RS iteration 1: Global Deviance = 502.7221 
GAMLSS-RS iteration 2: Global Deviance = 502.722 
GAMLSS-RS iteration 1: Global Deviance = 501.8834 
GAMLSS-RS iteration 2: Global Deviance = 501.8833 
GAMLSS-RS iteration 1: Global Deviance = 501.07 
GAMLSS-RS iteration 2: Global Deviance = 501.0699 
GAMLSS-RS iteration 1: Global Deviance = 500.2836 
GAMLSS-RS iteration 2: Global Deviance = 500.2834 
GAMLSS-RS iteration 1: Global Deviance = 499.526 
GAMLSS-RS iteration 2: Global Deviance = 499.5259 
GAMLSS-RS iteration 1: Global Deviance = 498.7993 
GAMLSS-RS iteration 2: Global Deviance = 498.7992 
GAMLSS-RS iteration 1: Global Deviance = 498.1053 
GAMLSS-RS iteration 2: Global Deviance = 498.1052 
GAMLSS-RS iteration 1: Global Deviance = 497.4458 
GAMLSS-RS iteration 2: Global Deviance = 497.4458 
GAMLSS-RS iteration 1: Global Deviance = 496.8227 
GAMLSS-RS iteration 2: Global Deviance = 496.8227 
GAMLSS-RS iteration 1: Global Deviance = 496.2377 
GAMLSS-RS iteration 2: Global Deviance = 496.2377 
GAMLSS-RS iteration 1: Global Deviance = 495.6924 
GAMLSS-RS iteration 2: Global Deviance = 495.6924 
GAMLSS-RS iteration 1: Global Deviance = 495.1883 
GAMLSS-RS iteration 2: Global Deviance = 495.1883 
GAMLSS-RS iteration 1: Global Deviance = 494.7269 
GAMLSS-RS iteration 2: Global Deviance = 494.7269 
GAMLSS-RS iteration 1: Global Deviance = 494.3093 
GAMLSS-RS iteration 2: Global Deviance = 494.3093 
GAMLSS-RS iteration 1: Global Deviance = 493.9368 
GAMLSS-RS iteration 2: Global Deviance = 493.9368 
GAMLSS-RS iteration 1: Global Deviance = 493.6101 
GAMLSS-RS iteration 2: Global Deviance = 493.6101 
GAMLSS-RS iteration 1: Global Deviance = 493.33 
GAMLSS-RS iteration 2: Global Deviance = 493.33 
GAMLSS-RS iteration 1: Global Deviance = 493.0972 
GAMLSS-RS iteration 2: Global Deviance = 493.0972 
GAMLSS-RS iteration 1: Global Deviance = 492.9118 
GAMLSS-RS iteration 2: Global Deviance = 492.9118 
GAMLSS-RS iteration 1: Global Deviance = 492.7741 
GAMLSS-RS iteration 2: Global Deviance = 492.7741 
GAMLSS-RS iteration 1: Global Deviance = 492.6839 
GAMLSS-RS iteration 2: Global Deviance = 492.6839 
GAMLSS-RS iteration 1: Global Deviance = 492.6412 
GAMLSS-RS iteration 2: Global Deviance = 492.6412 
GAMLSS-RS iteration 1: Global Deviance = 492.6455 
GAMLSS-RS iteration 2: Global Deviance = 492.6455 
GAMLSS-RS iteration 1: Global Deviance = 492.6962 
GAMLSS-RS iteration 2: Global Deviance = 492.6962 
GAMLSS-RS iteration 1: Global Deviance = 492.7926 
GAMLSS-RS iteration 2: Global Deviance = 492.7926 
GAMLSS-RS iteration 1: Global Deviance = 492.9339 
GAMLSS-RS iteration 2: Global Deviance = 492.9339 
GAMLSS-RS iteration 1: Global Deviance = 493.119 
GAMLSS-RS iteration 2: Global Deviance = 493.119 
GAMLSS-RS iteration 1: Global Deviance = 493.3469 
GAMLSS-RS iteration 2: Global Deviance = 493.3468 
GAMLSS-RS iteration 1: Global Deviance = 493.6164 
GAMLSS-RS iteration 2: Global Deviance = 493.6163 
GAMLSS-RS iteration 1: Global Deviance = 493.9262 
GAMLSS-RS iteration 2: Global Deviance = 493.9261 
GAMLSS-RS iteration 1: Global Deviance = 494.275 
GAMLSS-RS iteration 2: Global Deviance = 494.2748 
GAMLSS-RS iteration 1: Global Deviance = 494.6614 
GAMLSS-RS iteration 2: Global Deviance = 494.6612 
GAMLSS-RS iteration 1: Global Deviance = 495.0839 
GAMLSS-RS iteration 2: Global Deviance = 495.0837 
GAMLSS-RS iteration 1: Global Deviance = 495.5412 
GAMLSS-RS iteration 2: Global Deviance = 495.541 
GAMLSS-RS iteration 1: Global Deviance = 496.0317 
GAMLSS-RS iteration 2: Global Deviance = 496.0315 
GAMLSS-RS iteration 1: Global Deviance = 496.554 
GAMLSS-RS iteration 2: Global Deviance = 496.5537 
GAMLSS-RS iteration 1: Global Deviance = 497.1066 
GAMLSS-RS iteration 2: Global Deviance = 497.1063 
GAMLSS-RS iteration 1: Global Deviance = 497.688 
GAMLSS-RS iteration 2: Global Deviance = 497.6876 
GAMLSS-RS iteration 1: Global Deviance = 498.2967 
GAMLSS-RS iteration 2: Global Deviance = 498.2963 
GAMLSS-RS iteration 1: Global Deviance = 498.9313 
GAMLSS-RS iteration 2: Global Deviance = 498.9309 
GAMLSS-RS iteration 1: Global Deviance = 499.5905 
GAMLSS-RS iteration 2: Global Deviance = 499.5901 
GAMLSS-RS iteration 1: Global Deviance = 500.2728 
GAMLSS-RS iteration 2: Global Deviance = 500.2724 
GAMLSS-RS iteration 1: Global Deviance = 500.9769 
GAMLSS-RS iteration 2: Global Deviance = 500.9764 
GAMLSS-RS iteration 1: Global Deviance = 501.7015 
GAMLSS-RS iteration 2: Global Deviance = 501.701 
GAMLSS-RS iteration 1: Global Deviance = 502.4454 
GAMLSS-RS iteration 2: Global Deviance = 502.4449 
******************************************************************* 
Best estimate of the fixed parameter is  0.087 
with a Global Deviance equal to  492.6412  at position  28 
A  95 % Confidence interval is: ( 0.07458119 , 0.1008640 ) 
******************************************************************* 
> mod1<-quote(gamlss(y ~ cs(x,df=this) + qrt, data = aids, family = NBI))
> prof.term(mod1, min=1, max=15, step=1, criterion="IC")
GAMLSS-RS iteration 1: Global Deviance = 419.4861 
GAMLSS-RS iteration 2: Global Deviance = 423.8037 
GAMLSS-RS iteration 3: Global Deviance = 425.0068 
GAMLSS-RS iteration 4: Global Deviance = 424.9905 
GAMLSS-RS iteration 5: Global Deviance = 424.9851 
GAMLSS-RS iteration 6: Global Deviance = 424.9848 
GAMLSS-RS iteration 1: Global Deviance = 388.6842 
GAMLSS-RS iteration 2: Global Deviance = 391.3464 
GAMLSS-RS iteration 3: Global Deviance = 391.3969 
GAMLSS-RS iteration 4: Global Deviance = 391.3964 
GAMLSS-RS iteration 1: Global Deviance = 379.5593 
GAMLSS-RS iteration 2: Global Deviance = 379.8544 
GAMLSS-RS iteration 3: Global Deviance = 379.8779 
GAMLSS-RS iteration 4: Global Deviance = 379.8779 
GAMLSS-RS iteration 1: Global Deviance = 373.89 
GAMLSS-RS iteration 2: Global Deviance = 373.9232 
GAMLSS-RS iteration 3: Global Deviance = 373.9311 
GAMLSS-RS iteration 4: Global Deviance = 373.9313 
GAMLSS-RS iteration 1: Global Deviance = 369.3946 
GAMLSS-RS iteration 2: Global Deviance = 369.4118 
GAMLSS-RS iteration 3: Global Deviance = 369.4167 
GAMLSS-RS iteration 4: Global Deviance = 369.4169 
GAMLSS-RS iteration 1: Global Deviance = 365.4751 
GAMLSS-RS iteration 2: Global Deviance = 365.5071 
GAMLSS-RS iteration 3: Global Deviance = 365.5113 
GAMLSS-RS iteration 4: Global Deviance = 365.5116 
GAMLSS-RS iteration 1: Global Deviance = 362.0692 
GAMLSS-RS iteration 2: Global Deviance = 362.109 
GAMLSS-RS iteration 3: Global Deviance = 362.1121 
GAMLSS-RS iteration 4: Global Deviance = 362.1124 
GAMLSS-RS iteration 1: Global Deviance = 359.1964 
GAMLSS-RS iteration 2: Global Deviance = 359.2321 
GAMLSS-RS iteration 3: Global Deviance = 359.2344 
GAMLSS-RS iteration 4: Global Deviance = 359.2348 
GAMLSS-RS iteration 1: Global Deviance = 356.8284 
GAMLSS-RS iteration 2: Global Deviance = 356.8543 
GAMLSS-RS iteration 3: Global Deviance = 356.8559 
GAMLSS-RS iteration 4: Global Deviance = 356.8563 
GAMLSS-RS iteration 1: Global Deviance = 354.8884 
GAMLSS-RS iteration 2: Global Deviance = 354.9056 
GAMLSS-RS iteration 3: Global Deviance = 354.9066 
GAMLSS-RS iteration 4: Global Deviance = 354.9069 
GAMLSS-RS iteration 1: Global Deviance = 353.2625 
GAMLSS-RS iteration 2: Global Deviance = 353.2737 
GAMLSS-RS iteration 3: Global Deviance = 353.2744 
GAMLSS-RS iteration 1: Global Deviance = 351.8566 
GAMLSS-RS iteration 2: Global Deviance = 351.8638 
GAMLSS-RS iteration 3: Global Deviance = 351.8643 
GAMLSS-RS iteration 1: Global Deviance = 350.5468 
GAMLSS-RS iteration 2: Global Deviance = 350.551 
GAMLSS-RS iteration 3: Global Deviance = 350.5514 
GAMLSS-RS iteration 1: Global Deviance = 349.2648 
GAMLSS-RS iteration 2: Global Deviance = 349.2661 
GAMLSS-RS iteration 3: Global Deviance = 349.2667 
GAMLSS-RS iteration 1: Global Deviance = 347.9362 
GAMLSS-RS iteration 2: Global Deviance = 347.9333 
GAMLSS-RS iteration 3: Global Deviance = 347.9341 
> mod2 <- quote(gamlss(y ~ x+I((x>this)*(x-this))+qrt,family=NBI,data=aids))
> prof.term(mod2, min=1, max=45, step=1, criterion="GD")
GAMLSS-RS iteration 1: Global Deviance = 492.7247 
GAMLSS-RS iteration 2: Global Deviance = 492.6375 
GAMLSS-RS iteration 3: Global Deviance = 492.6373 
GAMLSS-RS iteration 1: Global Deviance = 483.1175 
GAMLSS-RS iteration 2: Global Deviance = 483.1128 
GAMLSS-RS iteration 3: Global Deviance = 483.1127 
GAMLSS-RS iteration 1: Global Deviance = 478.249 
GAMLSS-RS iteration 2: Global Deviance = 478.2475 
GAMLSS-RS iteration 3: Global Deviance = 478.2475 
GAMLSS-RS iteration 1: Global Deviance = 474.0472 
GAMLSS-RS iteration 2: Global Deviance = 474.0458 
GAMLSS-RS iteration 3: Global Deviance = 474.0458 
GAMLSS-RS iteration 1: Global Deviance = 469.3261 
GAMLSS-RS iteration 2: Global Deviance = 469.3245 
GAMLSS-RS iteration 3: Global Deviance = 469.3245 
GAMLSS-RS iteration 1: Global Deviance = 463.9904 
GAMLSS-RS iteration 2: Global Deviance = 463.9877 
GAMLSS-RS iteration 3: Global Deviance = 463.9877 
GAMLSS-RS iteration 1: Global Deviance = 457.2531 
GAMLSS-RS iteration 2: Global Deviance = 457.2481 
GAMLSS-RS iteration 3: Global Deviance = 457.2481 
GAMLSS-RS iteration 1: Global Deviance = 450.9153 
GAMLSS-RS iteration 2: Global Deviance = 450.9104 
GAMLSS-RS iteration 3: Global Deviance = 450.911 
GAMLSS-RS iteration 1: Global Deviance = 445.2555 
GAMLSS-RS iteration 2: Global Deviance = 445.2512 
GAMLSS-RS iteration 3: Global Deviance = 445.2512 
GAMLSS-RS iteration 1: Global Deviance = 439.8397 
GAMLSS-RS iteration 2: Global Deviance = 439.8355 
GAMLSS-RS iteration 3: Global Deviance = 439.8354 
GAMLSS-RS iteration 1: Global Deviance = 434.1987 
GAMLSS-RS iteration 2: Global Deviance = 434.1951 
GAMLSS-RS iteration 3: Global Deviance = 434.1951 
GAMLSS-RS iteration 1: Global Deviance = 427.9561 
GAMLSS-RS iteration 2: Global Deviance = 427.9518 
GAMLSS-RS iteration 3: Global Deviance = 427.9517 
GAMLSS-RS iteration 1: Global Deviance = 421.3125 
GAMLSS-RS iteration 2: Global Deviance = 421.3074 
GAMLSS-RS iteration 3: Global Deviance = 421.3074 
GAMLSS-RS iteration 1: Global Deviance = 412.7375 
GAMLSS-RS iteration 2: Global Deviance = 412.7268 
GAMLSS-RS iteration 3: Global Deviance = 412.7268 
GAMLSS-RS iteration 1: Global Deviance = 402.9311 
GAMLSS-RS iteration 2: Global Deviance = 402.9131 
GAMLSS-RS iteration 3: Global Deviance = 402.9129 
GAMLSS-RS iteration 1: Global Deviance = 392.778 
GAMLSS-RS iteration 2: Global Deviance = 392.7554 
GAMLSS-RS iteration 3: Global Deviance = 392.7553 
GAMLSS-RS iteration 1: Global Deviance = 382.8674 
GAMLSS-RS iteration 2: Global Deviance = 382.848 
GAMLSS-RS iteration 3: Global Deviance = 382.848 
GAMLSS-RS iteration 1: Global Deviance = 377.8746 
GAMLSS-RS iteration 2: Global Deviance = 377.8707 
GAMLSS-RS iteration 3: Global Deviance = 377.8706 
GAMLSS-RS iteration 1: Global Deviance = 378.9361 
GAMLSS-RS iteration 2: Global Deviance = 378.936 
GAMLSS-RS iteration 1: Global Deviance = 385.9809 
GAMLSS-RS iteration 2: Global Deviance = 385.9718 
GAMLSS-RS iteration 3: Global Deviance = 385.9718 
GAMLSS-RS iteration 1: Global Deviance = 394.6293 
GAMLSS-RS iteration 2: Global Deviance = 394.6019 
GAMLSS-RS iteration 3: Global Deviance = 394.6017 
GAMLSS-RS iteration 1: Global Deviance = 404.4438 
GAMLSS-RS iteration 2: Global Deviance = 404.3851 
GAMLSS-RS iteration 3: Global Deviance = 404.3848 
GAMLSS-RS iteration 1: Global Deviance = 413.5801 
GAMLSS-RS iteration 2: Global Deviance = 413.514 
GAMLSS-RS iteration 3: Global Deviance = 413.5135 
GAMLSS-RS iteration 1: Global Deviance = 421.1914 
GAMLSS-RS iteration 2: Global Deviance = 421.1412 
GAMLSS-RS iteration 3: Global Deviance = 421.1408 
GAMLSS-RS iteration 1: Global Deviance = 428.116 
GAMLSS-RS iteration 2: Global Deviance = 428.0718 
GAMLSS-RS iteration 3: Global Deviance = 428.0714 
GAMLSS-RS iteration 1: Global Deviance = 434.0359 
GAMLSS-RS iteration 2: Global Deviance = 434.0053 
GAMLSS-RS iteration 3: Global Deviance = 434.005 
GAMLSS-RS iteration 1: Global Deviance = 439.0664 
GAMLSS-RS iteration 2: Global Deviance = 439.0459 
GAMLSS-RS iteration 3: Global Deviance = 439.0456 
GAMLSS-RS iteration 1: Global Deviance = 443.8416 
GAMLSS-RS iteration 2: Global Deviance = 443.8229 
GAMLSS-RS iteration 3: Global Deviance = 443.8228 
GAMLSS-RS iteration 1: Global Deviance = 447.5496 
GAMLSS-RS iteration 2: Global Deviance = 447.5401 
GAMLSS-RS iteration 3: Global Deviance = 447.54 
GAMLSS-RS iteration 1: Global Deviance = 451.4804 
GAMLSS-RS iteration 2: Global Deviance = 451.4696 
GAMLSS-RS iteration 3: Global Deviance = 451.4695 
GAMLSS-RS iteration 1: Global Deviance = 455.4097 
GAMLSS-RS iteration 2: Global Deviance = 455.3994 
GAMLSS-RS iteration 3: Global Deviance = 455.3994 
GAMLSS-RS iteration 1: Global Deviance = 459.2939 
GAMLSS-RS iteration 2: Global Deviance = 459.2846 
GAMLSS-RS iteration 3: Global Deviance = 459.2845 
GAMLSS-RS iteration 1: Global Deviance = 462.9457 
GAMLSS-RS iteration 2: Global Deviance = 462.9385 
GAMLSS-RS iteration 3: Global Deviance = 462.9384 
GAMLSS-RS iteration 1: Global Deviance = 466.451 
GAMLSS-RS iteration 2: Global Deviance = 466.4452 
GAMLSS-RS iteration 3: Global Deviance = 466.4452 
GAMLSS-RS iteration 1: Global Deviance = 469.7953 
GAMLSS-RS iteration 2: Global Deviance = 469.7907 
GAMLSS-RS iteration 3: Global Deviance = 469.7907 
GAMLSS-RS iteration 1: Global Deviance = 472.6051 
GAMLSS-RS iteration 2: Global Deviance = 472.6025 
GAMLSS-RS iteration 3: Global Deviance = 472.6025 
GAMLSS-RS iteration 1: Global Deviance = 475.4832 
GAMLSS-RS iteration 2: Global Deviance = 475.4808 
GAMLSS-RS iteration 3: Global Deviance = 475.4807 
GAMLSS-RS iteration 1: Global Deviance = 477.9931 
GAMLSS-RS iteration 2: Global Deviance = 477.9915 
GAMLSS-RS iteration 3: Global Deviance = 477.9915 
GAMLSS-RS iteration 1: Global Deviance = 480.5113 
GAMLSS-RS iteration 2: Global Deviance = 480.51 
GAMLSS-RS iteration 3: Global Deviance = 480.51 
GAMLSS-RS iteration 1: Global Deviance = 482.8512 
GAMLSS-RS iteration 2: Global Deviance = 482.8503 
GAMLSS-RS iteration 1: Global Deviance = 485.0938 
GAMLSS-RS iteration 2: Global Deviance = 485.093 
GAMLSS-RS iteration 1: Global Deviance = 487.2029 
GAMLSS-RS iteration 2: Global Deviance = 487.2023 
GAMLSS-RS iteration 1: Global Deviance = 489.2268 
GAMLSS-RS iteration 2: Global Deviance = 489.2263 
GAMLSS-RS iteration 1: Global Deviance = 490.8663 
GAMLSS-RS iteration 2: Global Deviance = 490.866 
GAMLSS-RS iteration 1: Global Deviance = 492.6376 
GAMLSS-RS iteration 2: Global Deviance = 492.6373 
******************************************************************* 
Best estimate of the fixed parameter is  18 
with a Global Deviance equal to  377.8706  at position  18 
A  95 % Confidence interval is: ( 17.22821 , 19.39456 ) 
******************************************************************* 
> rm(mod,mod1,mod2)
> 
> 
> 
> cleanEx(); ..nameEx <- "ps"
> 
> ### * ps
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ps
> ### Title: Specify a Penalised Spline Fit in a GAMLSS Formula
> ### Aliases: ps
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> attach(aids)
> # fitting a smoothing cubic spline with 7 degrees of freedom
> # plus the a quarterly  effect  
> aids1<-gamlss(y~ps(x,df=7)+qrt,data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 374.8622 
GAMLSS-RS iteration 2: Global Deviance = 374.8622 
> plot(x,y)
> lines(x,fitted(aids1))
> rm(aids1)
> detach(aids)
> 
> 
> 
> cleanEx(); ..nameEx <- "ra"
> 
> ### * ra
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ra
> ### Title: Specify Simple Random Effect In A GAMLSS Formula
> ### Aliases: ra
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> attach(aids)
> # fitting a loess curve with span=0.4 plus the a quarterly  effect 
> aids1<-gamlss(y~lo(x,span=0.4)+qrt,data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 399.7146 
GAMLSS-RS iteration 2: Global Deviance = 399.747 
GAMLSS-RS iteration 3: Global Deviance = 399.7476 
> # now we string the quarterly  effect using random 
> aids2<-gamlss(y~lo(x,span=0.4)+ra(qrt,df=2),data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 411.5657 
GAMLSS-RS iteration 2: Global Deviance = 411.6035 
GAMLSS-RS iteration 3: Global Deviance = 411.6047 
GAMLSS-RS iteration 4: Global Deviance = 411.6103 
GAMLSS-RS iteration 5: Global Deviance = 411.6164 
GAMLSS-RS iteration 6: Global Deviance = 411.6215 
GAMLSS-RS iteration 7: Global Deviance = 411.6266 
GAMLSS-RS iteration 8: Global Deviance = 411.6306 
GAMLSS-RS iteration 9: Global Deviance = 411.6347 
GAMLSS-RS iteration 10: Global Deviance = 411.6376 
GAMLSS-RS iteration 11: Global Deviance = 411.6401 
GAMLSS-RS iteration 12: Global Deviance = 411.6425 
GAMLSS-RS iteration 13: Global Deviance = 411.6444 
GAMLSS-RS iteration 14: Global Deviance = 411.6464 
GAMLSS-RS iteration 15: Global Deviance = 411.648 
GAMLSS-RS iteration 16: Global Deviance = 411.6492 
GAMLSS-RS iteration 17: Global Deviance = 411.6504 
GAMLSS-RS iteration 18: Global Deviance = 411.651 
> plot(x,y)
> lines(x,fitted(aids1),col="red")
> lines(x,fitted(aids2),col="purple")
> rm(aids1,aids2)
> detach(aids)
> 
> 
> 
> cleanEx(); ..nameEx <- "random"
> 
> ### * random
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: random
> ### Title: Specify a simple random effect in a GAMLSS Formula
> ### Aliases: random
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> attach(aids)
> # fitting a loess curve with span=0.4 plus the a quarterly  effect 
> aids1<-gamlss(y~lo(x,span=0.4)+qrt,data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 399.7146 
GAMLSS-RS iteration 2: Global Deviance = 399.747 
GAMLSS-RS iteration 3: Global Deviance = 399.7476 
> # now we string the quarterly  effect using random 
> aids2<-gamlss(y~lo(x,span=0.4)+random(qrt,df=2),data=aids,family=PO) # 
GAMLSS-RS iteration 1: Global Deviance = 411.5793 
GAMLSS-RS iteration 2: Global Deviance = 411.6402 
GAMLSS-RS iteration 3: Global Deviance = 411.6505 
GAMLSS-RS iteration 4: Global Deviance = 411.6544 
GAMLSS-RS iteration 5: Global Deviance = 411.6554 
> plot(x,y)
> lines(x,fitted(aids1),col="red")
> lines(x,fitted(aids2),col="purple")
> rm(aids1,aids2)
> detach(aids)
> 
> 
> 
> cleanEx(); ..nameEx <- "rc"
> 
> ### * rc
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rc
> ### Title: Specify Random Coefficients In A GAMLSS Formula
> ### Aliases: rc
> ### Keywords: regression
> 
> ### ** Examples
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "refit"
> 
> ### * refit
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: refit
> ### Title: Refit a GAMLSS model
> ### Aliases: refit
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> refit(h)
GAMLSS-RS iteration 3: Global Deviance = 416.8014 

Family:  c("PO", "Poisson") 
Fitting method: RS() 

Call:  gamlss(formula = y ~ poly(x, 3) + qrt, family = PO, data = aids,      start.from = h, control = itercontrol) 

Mu Coefficients:
(Intercept)  poly(x, 3)1  poly(x, 3)2  poly(x, 3)3         qrt2         qrt3  
   4.814473     8.256092    -3.353310     0.938816    -0.156381     0.009986  
       qrt4  
  -0.115705  

Degrees of Freedom for the fit: 7 Residual Deg. of Freedom   38 
Global Deviance:     416.801 
            AIC:     430.801 
            SBC:     443.448 
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "rent"
> 
> ### * rent
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rent
> ### Title: Rent data
> ### Aliases: rent
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(rent)
> attach(rent)
> plot(Fl,R)
> 
> 
> 
> cleanEx(); ..nameEx <- "residuals.gamlss"
> 
> ### * residuals.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: residuals.gamlss
> ### Title: Extract Residuals from GAMLSS model
> ### Aliases: residuals.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=NBI, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 383.4573 
GAMLSS-RS iteration 2: Global Deviance = 381.7155 
GAMLSS-RS iteration 3: Global Deviance = 381.7145 
GAMLSS-RS iteration 4: Global Deviance = 381.7145 
> plot(aids$x,resid(h))
> plot(aids$x,resid(h,"sigma") )
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "rqres.plot"
> 
> ### * rqres.plot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rqres.plot
> ### Title: Plotting Randomized Quantile Residuals
> ### Aliases: rqres.plot
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids) # fitting a model from a discrete distribution 
> h<-gamlss(y~cs(x,df=7)+qrt, family=NBI, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 365.8195 
GAMLSS-RS iteration 2: Global Deviance = 361.9924 
GAMLSS-RS iteration 3: Global Deviance = 362.1084 
GAMLSS-RS iteration 4: Global Deviance = 362.1114 
GAMLSS-RS iteration 5: Global Deviance = 362.1123 
> plot(h)
*******************************************************************
	 Summary of the Randomised Quantile Residuals
                           mean   =  -0.01086206 
                       variance   =  0.9525586 
               coef. of skewness  =  -0.5769322 
               coef. of kurtosis  =  3.25145 
Filliben correlation coefficient  =  0.9853362 
*******************************************************************
> # plot qq- plots from 6 realization of the randomized quantile residuals
> rqres.plot(h) 
> # a qq-plot from the medians from 40 realizations
> rqres.plot(h,howmany=40,all=FALSE) # 
> 
> 
> 
> cleanEx(); ..nameEx <- "stepAIC"
> 
> ### * stepAIC
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: stepAIC
> ### Title: Choose a model by GAIC in a Stepwise Algorithm
> ### Aliases: stepAIC stepGAIC
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(usair)
> # fitting all variables linearly 
> mod1<-gamlss(y~., data=usair, family=GA)
GAMLSS-RS iteration 1: Global Deviance = 303.1604 
GAMLSS-RS iteration 2: Global Deviance = 303.1602 
> # find the best subset for the mu
> mod2<-stepAIC(mod1)
Distribution parameter:  mu 
Start:  AIC= 319.16 
 y ~ x1 + x2 + x3 + x4 + x5 + x6 

       Df    AIC
- x6    1 317.16
<none>    319.16
- x5    1 320.57
- x3    1 321.39
- x4    1 324.08
- x2    1 326.92
- x1    1 327.58

Step:  AIC= 317.16 
 y ~ x1 + x2 + x3 + x4 + x5 

       Df    AIC
<none>    317.16
- x3    1 319.39
- x4    1 322.48
- x5    1 324.14
- x2    1 324.92
- x1    1 336.11
> mod2$anova
Stepwise Model Path 
Analysis of Deviance Table

Initial
mu
 Model:
y ~ x1 + x2 + x3 + x4 + x5 + x6

Final
mu
 Model:
y ~ x1 + x2 + x3 + x4 + x5


  Step Df    Deviance Resid. Df Resid. Dev      AIC
1                            33   303.1602 319.1602
2 - x6  1 0.001712207        34   303.1619 317.1619
> # find the best subset for sigma
> mod3<-stepAIC(mod2, what="sigma", scope=~x1+x2+x3+x4+x5+x6)
Distribution parameter:  sigma 
Start:  AIC= 317.16 
 ~1 

       Df    AIC
+ x5    1 315.58
+ x4    1 315.96
<none>    317.16
+ x3    1 317.39
+ x6    1 317.53
+ x1    1 317.79
+ x2    1 318.61

Step:  AIC= 315.58 
 ~x5 

Warning in RS() : Algorithm RS has not yet converged
       Df    AIC
+ x1    1 314.08
<none>    315.58
+ x3    1 315.63
+ x6    1 316.63
+ x2    1 317.07
- x5    1 317.16
+ x4    1 317.77

Step:  AIC= 314.08 
 ~x5 + x1 

       Df    AIC
<none>    314.08
+ x3    1 314.48
+ x6    1 314.80
+ x2    1 315.38
+ x4    1 315.51
- x1    1 315.58
- x5    1 317.79
> mod3$anova
Stepwise Model Path 
Analysis of Deviance Table

Initial
sigma
 Model:
~1

Final
sigma
 Model:
~x5 + x1


  Step Df Deviance Resid. Df Resid. Dev      AIC
1                         34   303.1619 317.1619
2 + x5  1 3.584164        33   299.5777 315.5777
3 + x1  1 3.493242        32   296.0845 314.0845
> # now use the stepGAIC function
> # creating a scope from the usair model frame 
> gs<-gamlss.scope(model.frame(y~x1+x2+x3+x4+x5+x6, data=usair))
> gs 
$x1
~1 + x1 + cs(x1)

$x2
~1 + x2 + cs(x2)

$x3
~1 + x3 + cs(x3)

$x4
~1 + x4 + cs(x4)

$x5
~1 + x5 + cs(x5)

$x6
~1 + x6 + cs(x6)

> mod4<-gamlss(y~1, data=usair, family=GA)
GAMLSS-RS iteration 1: Global Deviance = 349.7146 
GAMLSS-RS iteration 2: Global Deviance = 349.7146 
> mod5<-stepGAIC(mod4,gs)
Distribution parameter:  mu 
Start:  y ~ 1; AIC= 353.7146 
Trial:  y ~  x1 + 1 + 1 + 1 + 1 + 1; AIC= 338.0354 
Trial:  y ~  1 + x2 + 1 + 1 + 1 + 1; AIC= 343.0487 
Trial:  y ~  1 + 1 + x3 + 1 + 1 + 1; AIC= 349.2046 
Trial:  y ~  1 + 1 + 1 + x4 + 1 + 1; AIC= 355.0206 
Trial:  y ~  1 + 1 + 1 + 1 + x5 + 1; AIC= 355.4256 
Trial:  y ~  1 + 1 + 1 + 1 + 1 + x6; AIC= 343.9733 
Step :  y ~ x1 ; AIC= 338.0354 

Trial:  y ~  cs(x1) + 1 + 1 + 1 + 1 + 1; AIC= 337.3823 
Trial:  y ~  x1 + x2 + 1 + 1 + 1 + 1; AIC= 328.781 
Trial:  y ~  x1 + 1 + x3 + 1 + 1 + 1; AIC= 332.942 
Trial:  y ~  x1 + 1 + 1 + x4 + 1 + 1; AIC= 339.6497 
Trial:  y ~  x1 + 1 + 1 + 1 + x5 + 1; AIC= 335.1107 
Trial:  y ~  x1 + 1 + 1 + 1 + 1 + x6; AIC= 335.9902 
Step :  y ~ x1 + x2 ; AIC= 328.781 

Trial:  y ~  cs(x1) + x2 + 1 + 1 + 1 + 1; AIC= 328.8904 
Trial:  y ~  x1 + cs(x2) + 1 + 1 + 1 + 1; AIC= 333.3071 
Trial:  y ~  x1 + x2 + x3 + 1 + 1 + 1; AIC= 325.665 
Trial:  y ~  x1 + x2 + 1 + x4 + 1 + 1; AIC= 327.6493 
Trial:  y ~  x1 + x2 + 1 + 1 + x5 + 1; AIC= 324.4414 
Trial:  y ~  x1 + x2 + 1 + 1 + 1 + x6; AIC= 325.5355 
Step :  y ~ x1 + x2 + x5 ; AIC= 324.4414 

Trial:  y ~  1 + x2 + 1 + 1 + x5 + 1; AIC= 344.2785 
Trial:  y ~  cs(x1) + x2 + 1 + 1 + x5 + 1; AIC= 323.8597 
Trial:  y ~  x1 + cs(x2) + 1 + 1 + x5 + 1; AIC= 328.7894 
Trial:  y ~  x1 + x2 + x3 + 1 + x5 + 1; AIC= 322.4756 
Trial:  y ~  x1 + x2 + 1 + x4 + x5 + 1; AIC= 319.3901 
Trial:  y ~  x1 + x2 + 1 + 1 + cs(x5) + 1; AIC= 316.6342 
Trial:  y ~  x1 + x2 + 1 + 1 + x5 + x6; AIC= 326.1607 
Step :  y ~ x1 + x2 + cs(x5) ; AIC= 316.6342 

Trial:  y ~  1 + x2 + 1 + 1 + cs(x5) + 1; AIC= 331.4525 
Trial:  y ~  cs(x1) + x2 + 1 + 1 + cs(x5) + 1; AIC= 318.5661 
Trial:  y ~  x1 + 1 + 1 + 1 + cs(x5) + 1; AIC= 331.739 
Trial:  y ~  x1 + cs(x2) + 1 + 1 + cs(x5) + 1; AIC= 319.0208 
Trial:  y ~  x1 + x2 + x3 + 1 + cs(x5) + 1; AIC= 315.4548 
Trial:  y ~  x1 + x2 + 1 + x4 + cs(x5) + 1; AIC= 316.2664 
Trial:  y ~  x1 + x2 + 1 + 1 + cs(x5) + x6; AIC= 318.2662 
Step :  y ~ x1 + x2 + x3 + cs(x5) ; AIC= 315.4548 

Trial:  y ~  1 + x2 + x3 + 1 + cs(x5) + 1; AIC= 323.81 
Trial:  y ~  cs(x1) + x2 + x3 + 1 + cs(x5) + 1; AIC= 317.3773 
Trial:  y ~  x1 + 1 + x3 + 1 + cs(x5) + 1; AIC= 321.4542 
Trial:  y ~  x1 + cs(x2) + x3 + 1 + cs(x5) + 1; AIC= 318.197 
Trial:  y ~  x1 + x2 + cs(x3) + 1 + cs(x5) + 1; AIC= 319.0917 
Trial:  y ~  x1 + x2 + x3 + x4 + cs(x5) + 1; AIC= 314.6209 
Trial:  y ~  x1 + x2 + x3 + 1 + cs(x5) + x6; AIC= 316.7411 
Step :  y ~ x1 + x2 + x3 + x4 + cs(x5) ; AIC= 314.6209 

Trial:  y ~  1 + x2 + x3 + x4 + cs(x5) + 1; AIC= 324.2394 
Trial:  y ~  cs(x1) + x2 + x3 + x4 + cs(x5) + 1; AIC= 316.4484 
Trial:  y ~  x1 + 1 + x3 + x4 + cs(x5) + 1; AIC= 321.8246 
Trial:  y ~  x1 + cs(x2) + x3 + x4 + cs(x5) + 1; AIC= 317.0621 
Trial:  y ~  x1 + x2 + cs(x3) + x4 + cs(x5) + 1; AIC= 318.0886 
Trial:  y ~  x1 + x2 + x3 + cs(x4) + cs(x5) + 1; AIC= 303.8334 
Trial:  y ~  x1 + x2 + x3 + x4 + x5 + 1; AIC= 317.1619 
Trial:  y ~  x1 + x2 + x3 + x4 + cs(x5) + x6; AIC= 316.4185 
Step :  y ~ x1 + x2 + x3 + cs(x4) + cs(x5) ; AIC= 303.8334 

Trial:  y ~  1 + x2 + x3 + cs(x4) + cs(x5) + 1; AIC= 311.8573 
Trial:  y ~  cs(x1) + x2 + x3 + cs(x4) + cs(x5) + 1; AIC= 305.8555 
Trial:  y ~  x1 + 1 + x3 + cs(x4) + cs(x5) + 1; AIC= 315.1206 
Trial:  y ~  x1 + cs(x2) + x3 + cs(x4) + cs(x5) + 1; AIC= 305.2431 
Trial:  y ~  x1 + x2 + 1 + cs(x4) + cs(x5) + 1; AIC= 309.6543 
Trial:  y ~  x1 + x2 + cs(x3) + cs(x4) + cs(x5) + 1; AIC= 307.9505 
Trial:  y ~  x1 + x2 + x3 + cs(x4) + x5 + 1; AIC= 305.1264 
Trial:  y ~  x1 + x2 + x3 + cs(x4) + cs(x5) + x6; AIC= 305.3408 
> mod5$anova
  From     To        Df   Deviance Resid. Df Resid. Dev      AIC
1                    NA         NA  39.00000   349.7146 353.7146
2          x1 -1.000000 -17.679187  38.00000   332.0354 338.0354
3          x2 -1.000000 -11.254434  37.00000   320.7810 328.7810
4          x5 -1.000000  -6.339651  36.00000   314.4414 324.4414
5   x5 cs(x5) -2.999638 -13.806449  33.00036   300.6349 316.6342
6          x3 -1.000000  -3.179355  32.00036   297.4555 315.4548
7          x4 -1.000000  -2.833942  31.00036   294.6216 314.6209
8   x4 cs(x4) -2.999290 -16.786090  28.00107   277.8355 303.8334
> mod6<-stepAIC(mod5, what="sigma", scope=~x1+x2+x3+x4+x5+x6)
Distribution parameter:  sigma 
Start:  AIC= 303.83 
 ~1 

Warning in RS() : Algorithm RS has not yet converged
            Df    AIC
+ x3   1.00198 299.93
+ x2   1.00098 302.09
<none>         303.83
+ x4   0.99944 304.37
+ x5   1.00240 305.59
+ x1   0.99962 305.63
+ x6   1.00108 305.83
Warning in RS() : Algorithm RS has not yet converged

Step:  AIC= 299.93 
 ~x3 

Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
            Df    AIC
+ x4   0.99774 299.50
<none>         299.93
+ x6   1.00029 300.85
+ x5   0.99906 301.22
+ x1   0.99731 302.09
+ x2   0.99900 302.79
- x3   1.00198 303.83
Warning in RS() : Algorithm RS has not yet converged

Step:  AIC= 299.5 
 ~x3 + x4 

Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
Warning in RS() : Algorithm RS has not yet converged
            Df    AIC
<none>         299.50
- x4   0.99774 299.93
+ x5   1.00124 300.12
+ x6   1.00211 301.17
+ x1   1.00112 301.25
+ x2   1.00243 302.38
- x3   1.00027 304.37
> mod6$anova
Stepwise Model Path 
Analysis of Deviance Table

Initial
sigma
 Model:
~1

Final
sigma
 Model:
~x3 + x4


  Step        Df Deviance Resid. Df Resid. Dev      AIC
1                          28.00107   277.8355 303.8334
2 + x3 1.0019760 5.908179  26.99910   271.9273 299.9291
3 + x4 0.9977373 2.427305  26.00136   269.5000 299.4973
> mod6

Family:  c("GA", "Gamma") 
Fitting method: RS() 

Call:  gamlss(formula = y ~ x1 + x2 + x3 + cs(x4) + cs(x5), sigma.formula = ~x3 +      x4, family = GA, data = usair, trace = FALSE) 

Mu Coefficients:
(Intercept)           x1           x2           x3       cs(x4)       cs(x5)  
  6.9192099   -0.0544108    0.0009063   -0.0004933   -0.1549188    0.0172195  
Sigma Coefficients:
(Intercept)           x3           x4  
  -2.288122    -0.001156     0.178630  

Degrees of Freedom for the fit: 14.99864 Residual Deg. of Freedom   26.00136 
Global Deviance:     269.5 
            AIC:     299.497 
            SBC:     325.199 
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "summary.gamlss"
> 
> ### * summary.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: summary.gamlss
> ### Title: Summarizes a GAMLSS fitted model
> ### Aliases: summary.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> h<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids) # 
GAMLSS-RS iteration 1: Global Deviance = 416.8014 
GAMLSS-RS iteration 2: Global Deviance = 416.8014 
> summary(h)
*******************************************************************
Family:  c("PO", "Poisson") 
Call:  gamlss(formula = y ~ poly(x, 3) + qrt, family = PO, data = aids) 
Fitting method: RS() 

-------------------------------------------------------------------
Mu link function:  log
Mu Coefficients:
              Estimate  Std. Error   z value   Pr(>|z|)
(Intercept)   4.814473     0.02821  170.6853  0.000e+00
poly(x, 3)1   8.256092     0.19436   42.4776  0.000e+00
poly(x, 3)2  -3.353310     0.16740  -20.0318  2.911e-89
poly(x, 3)3   0.938816     0.10937    8.5838  9.176e-18
qrt2         -0.156381     0.03060   -5.1110  3.205e-07
qrt3          0.009986     0.02868    0.3482  7.277e-01
qrt4         -0.115705     0.02911   -3.9745  7.053e-05

-------------------------------------------------------------------
No. of observations in the fit:  45 
Degrees of Freedom for the fit:  7
      Residual Deg. of Freedom:  38 
                      at cycle:  2 
 
Global Deviance:     416.8014 
            AIC:     430.8014 
            SBC:     443.4481 
*******************************************************************
> rm(h)
> 
> 
> 
> cleanEx(); ..nameEx <- "term.plot"
> 
> ### * term.plot
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: term.plot
> ### Title: Plot regression terms for a specified parameter of a GAMLSS
> ###   object
> ### Aliases: term.plot
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> a<-gamlss(y~cs(x,df=8)+qrt,data=aids,family=NBI)
GAMLSS-RS iteration 1: Global Deviance = 362.943 
GAMLSS-RS iteration 2: Global Deviance = 359.1257 
GAMLSS-RS iteration 3: Global Deviance = 359.229 
GAMLSS-RS iteration 4: Global Deviance = 359.2342 
GAMLSS-RS iteration 5: Global Deviance = 359.2348 
> term.plot(a, se=TRUE)
> rm(a)
> 
> 
> 
> cleanEx(); ..nameEx <- "update.gamlss"
> 
> ### * update.gamlss
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: update.gamlss
> ### Title: Update and Re-fit a GAMLSS Model
> ### Aliases: update.gamlss
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(aids)
> # fit a poisson model
> h.po <-gamlss(y~cs(x,2)+qrt, family=PO, data=aids) 
GAMLSS-RS iteration 1: Global Deviance = 448.1315 
GAMLSS-RS iteration 2: Global Deviance = 448.1315 
> # update with a negative binomial
> h.nb <-update(h.po, family=NBI) 
GAMLSS-RS iteration 1: Global Deviance = 388.5158 
GAMLSS-RS iteration 2: Global Deviance = 390.8177 
GAMLSS-RS iteration 3: Global Deviance = 391.3952 
GAMLSS-RS iteration 4: Global Deviance = 391.3994 
GAMLSS-RS iteration 5: Global Deviance = 391.3964 
GAMLSS-RS iteration 6: Global Deviance = 391.3962 
> # update the smoothing 
> h.nb1 <-update(h.nb,~cs(x,8)+qrt) 
GAMLSS-RS iteration 1: Global Deviance = 362.943 
GAMLSS-RS iteration 2: Global Deviance = 359.1257 
GAMLSS-RS iteration 3: Global Deviance = 359.229 
GAMLSS-RS iteration 4: Global Deviance = 359.2342 
GAMLSS-RS iteration 5: Global Deviance = 359.2348 
> # remove qrt
> h.nb2 <-update(h.nb1,~.-qrt)
GAMLSS-RS iteration 1: Global Deviance = 379.5511 
GAMLSS-RS iteration 2: Global Deviance = 379.5303 
GAMLSS-RS iteration 3: Global Deviance = 379.5628 
GAMLSS-RS iteration 4: Global Deviance = 379.5626 
> # put back qrt take log of y and fit a normal distribution 
> h.nb3 <-update(h.nb1,log(.)~.+qrt, family=NO) 
GAMLSS-RS iteration 1: Global Deviance = -42.3446 
GAMLSS-RS iteration 2: Global Deviance = -42.3446 
> # verify that it is the same 
> h.no<-gamlss(log(y)~cs(x,8)+qrt,data=aids )
GAMLSS-RS iteration 1: Global Deviance = -42.3446 
GAMLSS-RS iteration 2: Global Deviance = -42.3446 
> 
> 
> 
> cleanEx(); ..nameEx <- "usair"
> 
> ### * usair
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: usair
> ### Title: US air pollution data set
> ### Aliases: usair
> ### Keywords: datasets
> 
> ### ** Examples
> 
> data(usair)
> str(usair) 
`data.frame':	41 obs. of  7 variables:
 $ y : int  10 13 12 17 56 36 29 14 10 24 ...
 $ x1: num  70.3 61 56.7 51.9 49.1 54 57.3 68.4 75.5 61.5 ...
 $ x2: int  213 91 453 454 412 80 434 136 207 368 ...
 $ x3: int  582 132 716 515 158 80 757 529 335 497 ...
 $ x4: num  6 8.2 8.7 9 9 9 9.3 8.8 9 9.1 ...
 $ x5: num   7.05 48.52 20.66 12.95 43.37 ...
 $ x6: int  36 100 67 86 127 114 111 116 128 115 ...
> plot(usair)
> # a possible gamlss model
> ap<-gamlss(y~cs(x1,2)+x2+x3+cs(x4,2)+x5+cs(x6,3)+x4:x5, 
+                data=usair, family=GA(mu.link="inverse"))
GAMLSS-RS iteration 1: Global Deviance = 243.8559 
GAMLSS-RS iteration 2: Global Deviance = 243.8029 
GAMLSS-RS iteration 3: Global Deviance = 243.8024 
>  
> 
> 
> 
> cleanEx(); ..nameEx <- "wp"
> 
> ### * wp
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: wp
> ### Title: Worm plot
> ### Aliases: wp
> ### Keywords: regression
> 
> ### ** Examples
> 
> data(abdom)
> a<-gamlss(y~cs(x,df=3),sigma.fo=~cs(x,1),family=LO,data=abdom)
GAMLSS-RS iteration 1: Global Deviance = 4781.45 
GAMLSS-RS iteration 2: Global Deviance = 4781.096 
GAMLSS-RS iteration 3: Global Deviance = 4781.118 
GAMLSS-RS iteration 4: Global Deviance = 4781.118 
> wp(a)
> coeff1<-wp(a,abdom$x)
number of missing points from plot= 0  out of  154 
number of missing points from plot= 0  out of  153 
number of missing points from plot= 0  out of  156 
number of missing points from plot= 0  out of  147 
> coeff1
$classes
      [,1]  [,2]
[1,] 12.22 20.07
[2,] 20.07 27.07
[3,] 27.07 34.50
[4,] 34.50 42.50

$coef
            [,1]         [,2]        [,3]         [,4]
[1,]  0.03586757  0.043578914 0.019181573 -0.020052992
[2,] -0.09527334  0.008469508 0.053891087 -0.017925166
[3,]  0.03041428  0.058350522 0.030551341 -0.008392806
[4,] -0.02054224 -0.018977383 0.005199002  0.003507077

> rm(a,a1)
Warning: remove: variable "a1" was not found
> 
> 
> 
> ### * <FOOTER>
> ###
> cat("Time elapsed: ", proc.time() - get("ptime", env = .CheckExEnv),"\n")
Time elapsed:  404.44 6.53 417.3 0 0 
> grDevices::dev.off()
null device 
          1 
> ###
> ### Local variables: ***
> ### mode: outline-minor ***
> ### outline-regexp: "\\(> \\)?### [*]+" ***
> ### End: ***
> quit('no')