Which link function could be used for a glm where the response is per cent (0 - 100%)?

Question

I am thinking about building a model (glm) where the response variable (y) is the cover (in per cent) of a plant species in a defined area, dependant of environmental variables. However, I don't think it is a binomial link, as the output is not a probability between 0 and 1. What would you suggest to use as a link function?

Answer 1

Nice question!

You can start your modelling efforts by converting your cover variable so that it is expressed as a proportion rather than a percentage (i.e., simply divide the values of your original cover variable by 100).

Once you get the converted variable, look to see how many 0 and 1 values you have, if any, in addition to values falling in the (0,1) interval.

If you don't have any 0 or 1 values, you can use beta regression modelling. A beta regression model is not a glm model but it has similarities to a glm model. The beta regression model has a logit link function; this means that you are modelling the logit-transformed mean cover as a function of various predictor variables. (Recall that cover is now expressed as a proportion.)

If you have a substantial number of 0 values but no 1 values, you will need to use a zero-inflated beta regression model for your modelling.

If you have a substantial number of 1 values but no 0 values, you will need to use a one-inflated beta regression model.

If you have a substantial number of 0 values and a substantial number of 1 values, you will need to use a zero-and-one-inflated beta regression model.

Do you use R? If yes, the gamlss package is your best bet for implementing these kinds of models - see, for instance, page 107 of this document:

https://www.gamlss.com/wp-content/uploads/2013/01/book-2010-Athens1.pdf.

The inflated versions of the beta regression models have multiple links because they simultaneously model multiple parameters. Some of these links can be logit links, some can be log links, etc. In gamlss, the family distribution for these models can be one of:

1. BE for beta regression;

2. BEINF0 for zero-inflated beta regression;

3. BEINF1 for one-inflated beta regression;

4. BEINF for zero-and-one inflated beta regression.

Addendum

Both betareg() and gamlss() afford the ability to fit a beta regression model to your data and then use the model for prediction purposes. The R code below provides an example for how you can do this. Note that there are some differences in output between the model summaries produced by the two functions from the same dataset.

#===========================================================
# Load Gasoline data into R 
#===========================================================
library(betareg)
data("GasolineYield", package = "betareg")
#===========================================================
Fit beta regression model to Gasoline data using betareg()
#===========================================================
model_betareg <- betareg(yield ~ temp, 
                         link = "logit", 
                         link.phi = "log", 
                         data = GasolineYield)
summary(model_betareg)
#===========================================================
Fit beta regression model to Gasoline data using gamlss()
#===========================================================
library(gamlss)
model_gamlss <- gamlss(yield ~ temp, 
                       family = BE(mu.link = "logit", 
                                   sigma.link = "log"), 
                       data = GasolineYield)
summary(model_gamlss)
#===========================================================
Predict from model_betareg for a new temp value
#===========================================================
predict_betareg <- predict(model_betareg, 
                           newdata = data.frame(temp = 349), 
                           type="response")
predict_betareg
#===========================================================
Predict from model_gamlss for a new temp value
#===========================================================
predict_gamlss <- predict(model_gamlss, 
                          newdata = data.frame(temp = 349), 
                          type="response")
predict_gamlss

The yield value predicted by model_betareg is 0.2046289 and the yield value produced by model_gamlss is 0.204634. These are very close, as expected.

Note that the prediction is performed on the scale of the response variable yield, so that the predicted values are expressed as proportions in the interval (0,1).

The newdata must be a dataframe which lists the values of all of the predictor variables included in the right-hand side of the beta regression model. The type of these values must match the type of the predictors variables as seen by R via the str() command:

str(Gasoline)

Here's the output for this str() command:

str(Gasoline)
Classes ‘nfnGroupedData’, ‘nfGroupedData’, ‘groupedData’ and 'data.frame':      
32 obs. of  6 variables:
$ yield   : num  6.9 14.4 7.4 8.5 8 2.8 5 12.2 10 15.2 ...
$ endpoint: num  235 307 212 365 218 235 285 205 267 300 ...
$ Sample  : Ord.factor w/ 10 levels "1"<"2"<"7"<"9"<..: 8 9 5 1 3 4 7 10 6 8 
...
$ API     : num  38.4 40.3 40 31.8 40.8 41.3 38.1 50.8 32.2 38.4 ...
$ vapor   : num  6.1 4.8 6.1 0.2 3.5 1.8 1.2 8.6 5.2 6.1 ...
$ ASTM    : num  220 231 217 316 210 267 274 190 236 220 ...

If temp were listed as a factor in your dataset (e.g., "High", "Medium", "Low"), you would specify it as such in your newdata:

predict_gamlss <- predict(model_gamlss, 
                          newdata = data.frame(temp = c("High","Low")), 
                          type="response")
predict_gamlss