0
$\begingroup$

I want to identify if the response variable on the example data-set below is different between A and B groups.

DATA.SET
   GROUP   RESPONSE
     A    0.04965318
     A    0.02699422
     A    0.01849711
     A    0.05248555
     A    0.04398844
     A          0.01
     B    0.11479769
     B    0.04398844
     B    0.01283237
     B    0.02132948
     B    0.04682081
     B    0.01283237

Previously data inspection indicates that this sample came from a Beta distribution. So, using R with GAMLSS package, I try to analyse the above data set with the following lines:

library(gamlss)

gamlss(DATA.SET$RESPONSE~DATA.SET$GROUP, family=BE)

and these error and warnings come out:

Error in if (!family$y.valid(y)) stop("response variable out of range") : 
  missing value where TRUE/FALSE needed

In addition: Warning messages:
1: In Ops.factor(y, 0) : ‘>’ not meaningful for factors
2: In Ops.factor(y, 1) : ‘<’ not meaningful for factors

For the error, I don't know what is wrong, since beta distribution ranges from 0 to 1, but excluding these extremes values. The example data-set is correct about this.

Concerning the warnings, as far I know this indicates that I'm dealing with dummy variables, that are represented in many situations as 0 and 1. So I changed group to 0 and 1:

  DATA.SET
   GROUP   RESPONSE
      0    0.04965318
      0    0.02699422
      0    0.01849711
      0    0.05248555
      0    0.04398844
      0    0.01000000
      1    0.11479769
      1    0.04398844
      1    0.01283237
      1    0.02132948
      1    0.04682081
      1    0.01283237

Then I applied the same command line used before and everything goes fine without errors or warnings:

gamlss(DATA.SET$RESPONSE~DATA.SET$GROUP, family=BE)
GAMLSS-RS iteration 1: Global Deviance = -41.4784 
GAMLSS-RS iteration 2: Global Deviance = -47.7151 
GAMLSS-RS iteration 3: Global Deviance = -53.5619 
GAMLSS-RS iteration 4: Global Deviance = -57.0877 
GAMLSS-RS iteration 5: Global Deviance = -57.9639 
GAMLSS-RS iteration 6: Global Deviance = -58.0438 
GAMLSS-RS iteration 7: Global Deviance = -58.0476 
GAMLSS-RS iteration 8: Global Deviance = -58.0478 

Family:  c("BE", "Beta") 
Fitting method: RS() 

Call:  gamlss(formula = DATA.SET$RESPONSE ~ DATA.SET$GROUP,      family = BE) 

Mu Coefficients:
   (Intercept)  DATA.SET$GROUP  
      -3.26171         0.05667  
Sigma Coefficients:
(Intercept)  
     -1.882  

 Degrees of Freedom for the fit: 3 Residual Deg. of Freedom   9 
Global Deviance:     -58.0478 
            AIC:     -52.0478 
            SBC:     -50.5931

I read that many packages automatically identify dummy variables and proceed the analysis. I don't find any source stating that gamlss function do this, or more detailed tutorials on similar situations, so I'm not confident with the approach I took.

so...my doubt is if the way I deal this problem is correct? there are any sort of obviously wrong issues overlooked?

Also, I would greatly appreciate any sort of references about this GAMLSS issue!

$\endgroup$
1
$\begingroup$

I have put the data into a text file and everything works well I am not sure what is the problem?

> da<- read.table("~/.../DataExample.r", header=T)
> str(da)
'data.frame':   12 obs. of  2 variables:
 $ GROUP   : Factor w/ 2 levels "A","B": 1 1 1 1 1 1 2 2 2 2 ...
 $ RESPONSE: num  0.0497 0.027 0.0185 0.0525 0.044 ...
> library(gamlss)
> m1 <- gamlss(RESPONSE~GROUP, data=da, family=BE)
GAMLSS-RS iteration 1: Global Deviance = -41.4784 
GAMLSS-RS iteration 2: Global Deviance = -47.7151 
GAMLSS-RS iteration 3: Global Deviance = -53.5619 
GAMLSS-RS iteration 4: Global Deviance = -57.0877 
GAMLSS-RS iteration 5: Global Deviance = -57.9639 
GAMLSS-RS iteration 6: Global Deviance = -58.0438 
GAMLSS-RS iteration 7: Global Deviance = -58.0476 
GAMLSS-RS iteration 8: Global Deviance = -58.0478 

If you want to check what GAMLSS is doing with factors see below

> model.matrix(~GROUP, data=da)
   (Intercept) GROUPB
1            1      0
2            1      0
3            1      0
4            1      0
5            1      0
6            1      0
7            1      1
8            1      1
9            1      1
10           1      1
11           1      1
12           1      1
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$GROUP
[1] "contr.treatment"
$\endgroup$
1
  • $\begingroup$ Thanks for your fast reply @Mikis Stasinopoulos! I tried your lines and they did the job giving the same result I got when changed manually A-B to 0-1. So my doubt is solved since both ways give the same results, I'm doing right, only on a different way. About the error and warning on my first try, now I think there is something related to the way I create my data frame in R. I made it binding two vectors, not using read.table as you, maybe something on this way to construct data frames is incompatible with GAMLSS. $\endgroup$ – Emanuel Jan 29 '18 at 18:14

Not the answer you're looking for? Browse other questions tagged or ask your own question.