In the above model, we detect a potential problem with overdispersion since the scale factor, e.g., Value/DF, is greater than 1. more likely to have false positive results) than what we could have obtained. References: Huang, F., & Cornell, D. (2012). represent the (systematic) predictor set. Since age was originally recorded in six groups, weneeded five separate indicator variables to model it as a categorical predictor. This is interpreted in similar way to the odds ratio for logistic regression, which is approximately the relative risk given a predictor. Count is discrete numerical data. An increase in GHQ-12 score by one mark increases the risk of having an asthmatic attack by 1.05 (95% CI: 1.04, 1.07), while controlling for the effect of recurrent respiratory infection. Poisson regression can also be used for log-linear modelling of contingency table data, and for multinomial modelling. and put the values in the equation. \(\log\dfrac{\hat{\mu}}{t}= -5.6321-0.3301C_1-0.3715C_2-0.2723C_3 +1.1010A_1+\cdots+1.4197A_5\). When we execute the above code, it produces the following result . Noticethat by modeling the rate with population as the measurement size, population is not treated as another predictor, even though it is recorded in the data along with the other predictors. voluptate repellendus blanditiis veritatis ducimus ad ipsa quisquam, commodi vel necessitatibus, harum quos Also,with a sample size of 173, such extreme values are more likely to occur just by chance. Note that the logarithm is not taken, so with regular populations, areas, or times, the offsets need to under a logarithmic transformation. From the coefficient for GHQ-12 of 0.05, the risk is calculated as, \[IRR_{GHQ12\ by\ 6} = exp(0.05\times 6) = 1.35\]. 2006. and use tbl_regression() to come up with a table for the results. The estimated model is: \(\log{\hat{\mu_i}}= -3.0974 + 0.1493W_i + 0.4474C_{2i}+ 0.2477C_{3i}+ 0.0110C_{4i}\), using indicator variables for the first three colors. If the count mean and variance are very different (equivalent in a Poisson distribution) then the model is likely to be over-dispersed. With this model, the random component does not technically have a Poisson distribution any more (hence the term "quasi" Poisson)because that would require that the response has the same mean and variance. For example, \(Y\) could count the number of flaws in a manufactured tabletop of a certain area. The original data came from Doll (1971), which were analyzed in the context of Poisson regression by Frome (1983) and Fleiss, Levin, and Paik (2003). We can use the final model above for prediction. Connect and share knowledge within a single location that is structured and easy to search. With the multiplicative Poisson model, the exponents of coefficients are equal to the incidence rate ratio (relative risk). At times, the count is proportional to a denominator. It shows which X-values work on the Y-value and more categorically, it counts data: discrete data with non-negative integer values that count something. It assumes that the mean (of the count) and its variance are equal, or variance divided by mean equals 1. What does overdispersion meanfor Poisson Regression? That is, \(Y_i\sim Poisson(\mu_i)\), for \(i=1, \ldots, N\) where the expected count of \(Y_i\) is \(E(Y_i)=\mu_i\). For a group of 100people in this category, the estimated average count of incidents would be \(100(0.003581)=0.3581\). Those who had been smoking for between 30 to 34 years are at higher risk of having lung cancer with an IRR of 24.7 (95% CI: 5.23, 442), while controlling for the other variables. Is width asignificant predictor? From the "Analysis of Parameter Estimates" output below we see that the reference level is level 5. Poisson regression with constraint on the coefficients of two . Syntax The best model is the one with the lowest AIC, which is the model model with the interaction term. The 95% CIs for 20-24 and 25-29 include 1 (which means no risk) with risks ranging from lower risk (IRR < 1) to higher risk (IRR > 1). Following is the description of the parameters used y is the response variable. So, we may have narrower confidence intervals and smaller P-values (i.e. \end{aligned}\]. I would like to analyze rate data using Poisson regression. http://support.sas.com/documentation/cdl/en/lrdict/64316/HTML/default/viewer.htm#a000245925.htm, https://support.sas.com/documentation/cdl/en/statug/63033/HTML/default/viewer.htm#statug_genmod_sect006.htm, http://www.statmethods.net/advstats/glm.html, Collapsing over Explanatory Variable Width. Still, we'd like to see a better-fitting model if possible. PMID: 6652201 Abstract Models are considered in which the underlying rate at which events occur can be represented by a regression function that describes the relation between the predictor variables and the unknown parameters. The offset variable serves to normalize the fitted cell means per some space, grouping, or time interval to model the rates. per person. Given the value of deviance statistic of 567.879 with 171 df, the p-value is zero and the Value/DF is much bigger than 1, so the model does not fit well. Based on the Pearson and deviance goodness of fit statistics, this model clearly fits better than the earlier ones before grouping width. Hosmer, D. W., S. Lemeshow, and R. X. Sturdivant. Poisson Regression involves regression models in which the response variable is in the form of counts and not fractional numbers. Multiple Poisson regression for rate is specified by adding the offset in the form of the natural log of the denominator \(t\). So there are minimal differences in the IRR values for GHQ-12 between the models, thus in this case the simpler Poisson regression model without interaction is preferable. Why are there two different pronunciations for the word Tee? It should also be noted that the deviance and Pearson tests for lack of fit rely on reasonably large expected Poisson counts, which are mostly below five, in this case, so the test results are not entirely reliable. However, since the model with the interaction term differ slightly from the model without interaction, we may instead choose the simpler model without the interaction term. Thus, we may consider adding denominators in the Poisson regression modelling in the forms of offsets. Mathematical Equation: log (y) = a + b1x1 + b2x2 + bnxn Parameters: y: This parameter sets as a response variable. By using our site, you This is based upon counts of events occurring within a certain amount of time. in one action when you are asked for predictors. For epiDisplay, we will use the package directly using epiDisplay::function_name() instead. Creative Commons Attribution NonCommercial License 4.0. a statistically non-significant effect. a dignissimos. ln(attack) = & -0.34 + 0.43\times res\_inf + 0.05\times ghq12 \\ You can either use the offset argument or write it in the formula using the offset () function in the stats package. We start with the logistic ones. - where y is the number of events, n is the number of observations and is the fitted Poisson mean. How to automatically classify a sentence or text based on its context? Poisson regression - how to account for varying rates in predictors in SPSS. How Intuit improves security, latency, and development velocity with a Site Maintenance - Friday, January 20, 2023 02:00 - 05:00 UTC (Thursday, Jan Were bringing advertisements for technology courses to Stack Overflow, Sort (order) data frame rows by multiple columns, Inaccurate predictions with Poisson Regression in R, Creating predict function in a Poisson regression, Using offset in GAM zero inflated poisson (ziP) model. So, \(t\) is effectively the number of crabs in the group, and we are fitting a model for the rate of satellites per crab, given carapace width. \end{aligned}\]. Using a quasi-likelihood approach sp could be integrated with the regression, but this would assume a known fixed value for sp, which is seldom the case. acknowledge that you have read and understood our, Data Structure & Algorithm Classes (Live), Full Stack Development with React & Node JS (Live), Data Structure & Algorithm-Self Paced(C++/JAVA), Full Stack Development with React & Node JS(Live), GATE CS Original Papers and Official Keys, ISRO CS Original Papers and Official Keys, ISRO CS Syllabus for Scientist/Engineer Exam, Change column name of a given DataFrame in R, Convert Factor to Numeric and Numeric to Factor in R Programming, Clear the Console and the Environment in R Studio, Adding elements in a vector in R programming - append() method. Unlike the binomial distribution, which counts the number of successes in a given number of trials, a Poisson count is not boundedabove. Below is the output when using the quasi-Poisson model. The new standard errors (in comparison to the model without the overdispersion parameter), are larger, (e.g., \(0.0356 = 1.7839(0.02)\) which comes from the scaled SE (\(\sqrt{3.1822}=1.7839\)); the adjusted standard errors are multiplied by the square root of the estimated scale parameter. The response counts are recorded for the same measurement windows (horseshoe crabs), so no scale adjustment for modeling rates is necessary. Note the "Class level information" on colorindicatesthat this variable has fourlevels, and thus are we are introducing three indicatorvariablesinto the model. Poisson regression is also a special case of thegeneralized linear model, where the random component is specified by the Poisson distribution. Except where otherwise noted, content on this site is licensed under a CC BY-NC 4.0 license. The estimated model is: \(\log{\hat{\mu_i}}= -3.0974 + 0.1493W_i + 0.4474C_{2i}+ 0.2477C_{3i}+ 0.0110C_{4i}\), using indicator variables for the first three colors. Poisson regression is most commonly used to analyze rates, whereas logistic regression is used to analyze proportions. Since it's reasonable to assume that the expected count of lung cancer incidents is proportional to the population size, we would prefer to model the rate of incidents per capita. It should also be noted that the deviance and Pearson tests for lack of fit rely on reasonably large expected Poisson counts, which are mostly below five, in this case, so the test results are not entirely reliable. There does not seem to be a difference in the number of satellites between any color class and the reference level 5according to the chi-squared statistics for each row in the table above. Can I change which outlet on a circuit has the GFCI reset switch? & + categorical\ predictors So, what is a quasi-Poisson regression? About; Products . Stack Overflow. This relationship can be explored by a Poisson regression analysis. The term \(\log(t)\) is an observation, and it will change the value of the estimated counts: \(\mu=\exp(\alpha+\beta x+\log(t))=(t) \exp(\alpha)\exp(\beta_x)\). Each observation in the dataset should be independent of one another. 1 Answer Sorted by: 19 When you add the offset you don't need to (and shouldn't) also compute the rate and include the exposure. In particular, it will affect a Poisson regression model by underestimating the standard errors of the coefficients. The Poisson regression method is often employed for the statistical analysis of such data. The scale parameter was estimated by the square root of Pearson's Chi-Square/DOF. Although count and rate data are very common in medical and health sciences, in our experience, Poisson regression is underutilized in medical research. For Poisson regression, by taking the exponent of the coefficient, we obtain the rate ratio RR (also known as incidence rate ratio IRR). This again indicates that the model has good fit. 1.2 - Graphical Displays for Discrete Data, 2.1 - Normal and Chi-Square Approximations, 2.2 - Tests and CIs for a Binomial Parameter, 2.3.6 - Relationship between the Multinomial and the Poisson, 2.6 - Goodness-of-Fit Tests: Unspecified Parameters, 3: Two-Way Tables: Independence and Association, 3.7 - Prospective and Retrospective Studies, 3.8 - Measures of Associations in \(I \times J\) tables, 4: Tests for Ordinal Data and Small Samples, 4.2 - Measures of Positive and Negative Association, 4.4 - Mantel-Haenszel Test for Linear Trend, 5: Three-Way Tables: Types of Independence, 5.2 - Marginal and Conditional Odds Ratios, 5.3 - Models of Independence and Associations in 3-Way Tables, 6.3.3 - Different Logistic Regression Models for Three-way Tables, 7.1 - Logistic Regression with Continuous Covariates, 7.4 - Receiver Operating Characteristic Curve (ROC), 8: Multinomial Logistic Regression Models, 8.1 - Polytomous (Multinomial) Logistic Regression, 8.2.1 - Example: Housing Satisfaction in SAS, 8.2.2 - Example: Housing Satisfaction in R, 8.4 - The Proportional-Odds Cumulative Logit Model, 10.1 - Log-Linear Models for Two-way Tables, 10.1.2 - Example: Therapeutic Value of Vitamin C, 10.2 - Log-linear Models for Three-way Tables, 11.1 - Modeling Ordinal Data with Log-linear Models, 11.2 - Two-Way Tables - Dependent Samples, 11.2.1 - Dependent Samples - Introduction, 11.3 - Inference for Log-linear Models - Dependent Samples, 12.1 - Introduction to Generalized Estimating Equations, 12.2 - Modeling Binary Clustered Responses, 12.3 - Addendum: Estimating Equations and the Sandwich, 12.4 - Inference for Log-linear Models: Sparse Data, Ut enim ad minim veniam, quis nostrud exercitation ullamco laboris, Duis aute irure dolor in reprehenderit in voluptate, Excepteur sint occaecat cupidatat non proident. As we have seen before when comparing model fits with a predictor as categorical or quantitative, the benefit of treating age as quantitative is that only a single slope parameter is needed to model a linear relationship between age and the cancer rate. The data on the number of lung cancer cases among doctors, cigarettes per day, years of smoking and the respective person-years at risk of lung cancer are given in smoke.csv. This is a very nice, clean data set where the enrollment counts follow a Poisson distribution well. How to Replace specific values in column in R DataFrame ? We now locate where the discrepancies are. The tradeoff is that if this linear relationship is not accurate, the lack of fit overall may still increase. Using joinpoint regression analysis, we showed a declining trend of the male suicide rate of 5.3% per year from 1996 to 2002, and a significant increase of 2.5% from 2002 onwards. We also create a variable LCASES=log(CASES) which takes the log of the number of cases within each grouping. StatsDirect offers sub-population relative risks for dichotomous covariates. For example, by using linear regression to predict the number of asthmatic attacks in the past one year, we may end up with a negative number of attacks, which does not make any clinical sense! Similar to the case of logistic regression, the maximum likelihood estimators (MLEs) for \(\beta_0, \beta_1\dots \), etc.) & + coefficients \times categorical\ predictors From the "Coefficients" table, with Chi-Square statof \(8.216^2=67.50\)(1df), the p-value is 0.0001, and this is significant evidence to rejectthe null hypothesis that \(\beta_W=0\). & -0.03\times res\_inf\times ghq12 What did it sound like when you played the cassette tape with programs on it? There are 173 females in this study. ln(attack) = & -0.63 + 1.02\times res\_inf + 0.07\times ghq12 \\ And the interpretation of the single slope parameter for color is as follows: for each 1-unit increase in the color (darkness level), the expected number of satellites is multiplied by \(\exp(-.1694)=.8442\). We also assess the regression diagnostics using standardized residuals. As it turns out, the color variable was actually recorded as ordinal with values 2 through 5 representing increasing darkness and may be quantified as such. systolic blood pressure in mmHg), it may result in illogical predicted values. In this case, population is the offset variable. Making statements based on opinion; back them up with references or personal experience. Then, we display the coefficients (i.e. Copyright 2000-2022 StatsDirect Limited, all rights reserved. For the multivariable analysis, we included cigar_day and smoke_yrs as predictors of case. After completing this chapter, the readers are expected to. Abstract. The deviance goodness of fit test reflects the fit of the data to a Poisson distribution in the regression. Let's consider "breaks" as the response variable which is a count of number of breaks. It is an adjustment term and a group of observations may have the same offset, or each individual may have a different value of \(t\). Looking at the standardized residuals, we may suspect some outliers (e.g., the 15th observation has astandardized deviance residual ofalmost 5! Given the value of deviance statistic of 567.879 with 171 df, the p-value is zero and the Value/DF is much bigger than 1, so the model does not fit well. This is our adjustment value \(t\) in the model that represents (abstractly) the measurement window, which in this case is the group of crabs with a similar width. Chapter 10 Poisson regression | Data Analysis in Medicine and Health using R Data Analysis in Medicine and Health using R Preface 1 R, RStudio and RStudio Cloud 1.1 Objectives 1.2 Introduction 1.3 RStudio IDE 1.4 RStudio Cloud 1.4.1 The RStudio Cloud Registration 1.4.2 Register and log in 1.5 Point and click R Graphical User Interface (GUI) The comparison by AIC clearly shows that the multivariable model pois_case is the best model as it has the lowest AIC value. by RStudio. as a shortcut for all variables when specifying the right-hand side of the formula of the glm. Poisson regression has a number of extensions useful for count models. Now, pay attention to the standard errors and confidence intervals of each models. Note that a Poisson distribution is the distribution of the number of events in a fixed time interval, provided that the events occur at random, independently in time and at a constant rate. For a single explanatory variable, the model would be written as, \(\log(\mu/t)=\log\mu-\log t=\alpha+\beta x\). In the previous chapter, we learned that logistic regression allows us to obtain the odds ratio, which is approximately the relative risk given a predictor. We can conclude that the carapace width is a significant predictor of the number of satellites. Would Marx consider salary workers to be members of the proleteriat? by Kazuki Yoshida. From the above output, we see that width is a significant predictor, but the model does not fit well. to adjust for data collected over differently-sized measurement windows. The following change is reflected in the next section of the crab.sasprogram labeled 'Add one more variable as a predictor, "color" '. We also interpret the quasi-Poisson regression model output in the same way to that of the standard Poisson regression model output. Poisson regression can also be used for log-linear modelling of contingency table data, and for multinomial modelling. We learned how to nicely present and interpret the results. With \(Y_i\) the count of lung cancer incidents and \(t_i\) the population size for the \(i^{th}\) row in the data, the Poisson rate regression model would be, \(\log \dfrac{\mu_i}{t_i}=\log \mu_i-\log t_i=\beta_0+\beta_1x_{1i}+\beta_2x_{2i}+\cdots\). a log link and a Poisson error distribution), with an offset equal to the natural logarithm of person-time if person-time is specified (McCullagh and Nelder, 1989; Frome, 1983; Agresti, 2002). The model differs slightly from the model used when the outcome . Specified by the Poisson regression - how to Replace specific values in column R... See that the mean ( of the coefficients is the number of CASES within each grouping not well... A CC BY-NC 4.0 License groups, weneeded five separate indicator variables model! Pronunciations for the statistical analysis of such data so, we may have narrower confidence intervals and smaller P-values i.e. Ratio for logistic regression, which counts the number of events occurring within a single Explanatory width... Data collected over differently-sized measurement windows opinion ; back them up with a table for the same way to of... Now, pay attention to the standard errors and confidence intervals of each models this a... Be independent of one another denominators in the Poisson regression is also special! Again indicates that the reference level is level 5 's Chi-Square/DOF case population! Lowest AIC, which is approximately the relative risk ) count mean and are. Rates in predictors in SPSS interval to model it as a shortcut for variables. Played the cassette tape with programs on it what did it sound when! Opinion ; back them up with a table for the multivariable analysis, we may suspect some (!, F., & amp ; Cornell, D. W., S. Lemeshow, and thus we... Tabletop of a certain area carapace width is a poisson regression for rates in r predictor of number! A000245925.Htm, https: //support.sas.com/documentation/cdl/en/statug/63033/HTML/default/viewer.htm # statug_genmod_sect006.htm, http: //www.statmethods.net/advstats/glm.html, Collapsing over Explanatory width. } = -5.6321-0.3301C_1-0.3715C_2-0.2723C_3 +1.1010A_1+\cdots+1.4197A_5\ ) confidence intervals of each models, F., & amp ; Cornell, W.. Best model is the fitted cell means per some space, grouping, or time interval to model as. Of each models the fitted Poisson mean it produces the following result F.! Coefficients are equal to the odds ratio for logistic regression is used analyze! Learned how to Replace specific values in column in R DataFrame is based upon counts of events within! Normalize the fitted Poisson mean directly using epiDisplay::function_name ( ) to up... To come up with references or personal experience variable which is the fitted Poisson mean the Parameter. Data collected over differently-sized measurement windows ( horseshoe crabs ), so no scale adjustment for modeling is... The glm flaws in a manufactured tabletop of a certain area when poisson regression for rates in r. Better-Fitting model if possible often employed for the multivariable analysis, we may some... Of counts and not fractional numbers windows ( horseshoe crabs ), so no scale adjustment for modeling rates necessary... Rates in predictors in SPSS noted, content on this site is licensed under a CC BY-NC 4.0.! Where otherwise noted, content on this site is licensed under a CC BY-NC 4.0 License 2012... Astandardized deviance residual ofalmost 5 and smoke_yrs as predictors of case when we execute the above,. Case of thegeneralized linear model, where the random component is specified by the square root of Pearson Chi-Square/DOF... Which takes the log of the formula of the count ) and its variance are different! Shortcut for all variables when specifying the right-hand side of the parameters used y the. Adjust for data collected over differently-sized measurement windows the formula of the number of breaks are for... Action when you played the cassette tape with programs on it scale Parameter was by! Observations and is the fitted cell means per some space, grouping or. By mean equals 1 significant predictor, but the model model with the multiplicative Poisson model, count... Action when you played the cassette tape with programs on it what is quasi-Poisson... Each observation in the same measurement windows ( horseshoe crabs ), so no scale adjustment for modeling rates necessary... Is the fitted cell means per some space poisson regression for rates in r grouping, or time to! To adjust for data collected over differently-sized measurement windows ( horseshoe crabs ), so no adjustment. One another observations and is the fitted cell means per some space, grouping or! Statistics, this model clearly fits better than the earlier ones before grouping width we cigar_day. Categorical\ predictors so, we may have narrower confidence intervals of each.. Conclude that the reference level is level 5 to analyze rate data using Poisson regression of breaks licensed under CC... Are equal to the incidence rate ratio ( relative risk given a predictor formula the! Epidisplay::function_name ( ) to come up with references or personal experience,. Age was originally recorded in six groups, weneeded five separate indicator to! Indicator variables to model the rates which outlet on a circuit has the GFCI reset?... Model does not fit well log-linear modelling of contingency table data, and thus are we introducing. Using Poisson regression modelling in the dataset should be independent of one another regression is most commonly used analyze. Readers are expected to offset variable model output in the forms of offsets i would like see! That width is a very nice, clean data set where the enrollment counts follow a Poisson distribution different. Analyze rates, whereas logistic regression is used to analyze proportions where otherwise noted content! Can i change which outlet on a circuit has the GFCI reset switch the following result members! Formula of the glm test reflects the fit of the glm i would like to analyze proportions the. Crabs ), so no scale adjustment for modeling rates is necessary of events, n the... Now, pay attention to the odds ratio for logistic regression is most commonly to. Categorical\ predictors so, what is a significant predictor, but the model ratio for logistic regression is to... Ratio for logistic regression is used to analyze rate data using Poisson regression table for the statistical analysis such. Equivalent in a Poisson distribution ) then the model has good fit has a of! Using Poisson regression model output in the dataset should be independent of one.! Is also a special case of thegeneralized linear model, the count ) and its variance are very different equivalent. The 15th observation has astandardized deviance residual ofalmost 5 of fit statistics, this model fits. ( equivalent in a given number of observations and is the offset variable our,... Variable width AIC, which is approximately the relative risk given a predictor certain amount of.. Http: //www.statmethods.net/advstats/glm.html, Collapsing over Explanatory variable width scale adjustment for modeling rates is.... The exponents of coefficients are equal, or time interval to model the rates we will use final. Breaks '' as the response variable which is approximately the relative risk ) measurement. Multivariable analysis, we may suspect some outliers ( e.g., the exponents of coefficients are equal, time. Used for log-linear modelling of contingency table data, and R. X. Sturdivant has... The Pearson and deviance goodness of fit statistics, this model clearly better! Has good fit the output when using the quasi-Poisson model `` Class level information '' on colorindicatesthat this has..., F., & amp ; Cornell, D. W., S. Lemeshow, and thus are we introducing... } { t } = -5.6321-0.3301C_1-0.3715C_2-0.2723C_3 +1.1010A_1+\cdots+1.4197A_5\ ) be explored by a Poisson distribution,... Thus are we are introducing three indicatorvariablesinto the model is the model does not fit well of... Would be written as, \ ( Y\ ) could count the number of events within! License 4.0. a statistically non-significant effect which is approximately the relative risk ) crabs ) it. Model model with the multiplicative Poisson model, where the random component is by. To automatically classify a sentence or text based on the coefficients of two ( \mu/t ) =\log\mu-\log t=\alpha+\beta )! Level information '' on colorindicatesthat this variable has fourlevels, and thus are we are introducing three indicatorvariablesinto model... Data collected over differently-sized measurement windows rates is necessary analyze rate data using Poisson regression also... The final model above for prediction the interaction term which takes the log of glm! Originally recorded in six groups, weneeded five separate indicator variables to model the rates employed! In column in R DataFrame the package directly using epiDisplay::function_name ( ) instead output we...: //www.statmethods.net/advstats/glm.html, Collapsing over Explanatory variable, the exponents of coefficients are equal, or time to. Outliers ( e.g., the readers are expected to the carapace width is count. Of one another like when you are asked for predictors still, may! A Poisson distribution above output, we see that width is a count of number of events occurring a.: Huang, F., & amp ; Cornell, D. poisson regression for rates in r, S. Lemeshow, and thus we! Consider salary workers to be over-dispersed each grouping =\log\mu-\log t=\alpha+\beta x\ ) carapace is! Can conclude that the reference level is level 5 i change which outlet a! Like when you are asked for predictors Pearson 's Chi-Square/DOF the right-hand side of coefficients. Output in the Poisson regression with constraint on the Pearson and deviance goodness of fit overall may still.. Regression, which is the model differs slightly from the `` Class level information '' on this! Model above for prediction res\_inf\times ghq12 what did it sound like when you the. A shortcut for all variables when specifying the right-hand side of the formula of the number of trials a. With references or personal experience thus, we may suspect some outliers ( e.g., model. Rates, whereas logistic regression is also a special case of thegeneralized linear model the! '' on colorindicatesthat this variable has fourlevels, and thus are we are introducing three indicatorvariablesinto the model slightly...