# Multiple regression with categorical and numeric predictors

I am relatively new to R, and I am trying to fit a model to data that consists of a categorical column and a numeric (integer) column. The dependent variable is a continuous number.

The data has the following format:

predCateg, predIntNum, ResponseVar

The data looks something like this:

ranking, age_in_years, wealth_indicator
category_A, 99, 1234.56
category_A, 21, 12.34
category_A, 42, 234.56
....
category_N, 105, 77.27


How would I model this (presumably, using a GLM), in R?

[]

It has just occurred to me (after analysing the data more thoroughly), that the categorical independent variable is in fact ordered. I have therefore modified the answer provided earlier as follows:

> fit2 <- glm(wealth_indicator ~ ordered(ranking) + age_in_years, data=amort2)
>
> fit2

Call:  glm(formula = wealth_indicator ~ ordered(ranking) + age_in_years,
data = amort2)

Coefficients:
(Intercept)  ordered(ranking).L  ordered(ranking).Q  ordered(ranking).C      age_in_years
0.0578500         -0.0055454         -0.0013000          0.0007603          0.0036818

Degrees of Freedom: 39 Total (i.e. Null);  35 Residual
Null Deviance:      0.004924
Residual Deviance: 0.00012      AIC: -383.2
>
> fit3 <- glm(wealth_indicator ~ ordered(ranking) + age_in_years + ordered(ranking)*age_in_years, data=amort2)
> fit3

Call:  glm(formula = wealth_indicator ~ ordered(ranking) + age_in_years +
ordered(ranking) * age_in_years, data = amort2)

Coefficients:
(Intercept)                ordered(ranking).L                ordered(ranking).Q
0.0578500                       -0.0018932                       -0.0039667
ordered(ranking).C                    age_in_years  ordered(ranking).L:age_in_years
0.0021019                        0.0036818                       -0.0006640
ordered(ranking).Q:age_in_years  ordered(ranking).C:age_in_years
0.0004848                       -0.0002439

Degrees of Freedom: 39 Total (i.e. Null);  32 Residual
Null Deviance:      0.004924
Residual Deviance: 5.931e-05    AIC: -405.4


I am a bit confused by what ordered(ranking).C, ordered(ranking).Q and ordered(ranking).L mean in the output, and would appreciate some help in understanding this output, and how to use it to predict the response variable.

Try this:

fit <- glm(wealth_indicator ~
factor(ranking) + age_in_years + factor(ranking) * age_in_years)


The factor() command will make sure that R knows that your variable is categorical. This is especially useful if your categories are indicated by integers, otherwise glm will interpret the variable as continuous.

The factor(ranking) * age_in_years term lets R know that you want to include the interaction term.

• why factor(ranking) and not as.factor(ranking) ? Mar 3, 2014 at 18:26
• I habitually use factor(x) so that I can include the levels argument if I wish. You could also use as.factor(x) if you wish, and it may in fact be faster, but I would think you'd need quite a large dataset for the speed of these functions to matter. Mar 3, 2014 at 18:35
• @PSchnell +1 for a concise answer and a good explanation. One thing however, is that I have realized that the factor is ordered. I modified the formula you provided slightly, to reflect the fact this (please see my edited question). However, I am not sure I understand the model output (specifically ordered(ranking).C, ordered(ranking).Q and ordered(ranking).L - what do they mean, and how do I use that to predict the response variable?) - any help will be much appreciated. Thanks Mar 4, 2014 at 5:36
• .L, .Q, and .C are, respectively, the coefficients for the ordered factor coded with linear, quadratic, and cubic contrasts. The command contr.poly(4) will show you the contrast matrix for an ordered factor with 4 levels (3 degrees of freedom, which is why you get up to a third order polynomial). contr.poly(4)[2, '.L'] will tell you what to plug in for the second ordered level in the linear term. Be aware that this assumes that it makes sense to consider the levels as equally spaced. If it doesn't, code your own contrast matrix. Mar 4, 2014 at 13:26
• According to ats.ucla.edu/stat/r/library/contrast_coding.htm#ORTHOGONAL, in R, when a continuous variable is modeled in response to an ordered factor variable, the default contrast applied is orthogonal polynomial coding onlinecourses.science.psu.edu/stat502/node/203. stats.stackexchange.com/a/206345/90600 has some rather extensive explanation. stats.stackexchange.com/a/207128/90600 goes even deeper. Apr 22, 2016 at 6:54

I recently answered Continuous dependent variable with ordinal independent variable, recommending the ordSmooth function in the ordPens package. This uses penalized regression to smooth dummy coefficients across levels of an ordinal variable so that they don't vary too greatly from one level to the next. E.g., you probably wouldn't want category_A's coefficient to be much more different from category_B than from category_N. You'd probably want coefficients to rise or fall monotonically, or at least change gradually across ranks. My answer to the linked question lists references for further info.

ordSmooth can also accommodate continuous (and nominal) variables; in your case, code could be:

SmoothFit=with(amort2,
ordSmooth(as.numeric(ordered(ranking)),wealth_indicator,z=age_in_years,lambda=.001))


ordSmooth requires numeric input for ordinal data, hence the as.numeric(ordered()) reformatting. z is for a numeric vector/matrix/data.frame of continuous predictors. lambda is the smoothing coefficient – larger values will push your coefficients closer to zero. Might be wise to start small here. Printing SmoothFit will give you coefficients and fitted values, but I'm afraid it leaves the rest to you.

In your method, the ordered(ranking).C/.Q/.L coefficients appear to be labeled as cubic, quadratic, and linear, respectively. If you try glm(rnorm(10)~ordered(rep(1:5,2))), you'll get a coefficient for ordered(rep(1:5, 2))^4 as well. I'm not really sure why these are denoted with exponents; I don't think it's modeling these as polynomial terms, because the coefficients are different for glm(y~x+I(x^2)+I(x^3)+I(x^4)) and scaled variants of this. They should be basic dummy codes.