# Consequence of choosing wrong functional of covariates in GLM/GAM

I'm modelling the mood of teenagers in a really big school. Response is 'good mood' and 'bad mood'. One of the variables that is used to explain the students mood is "Area of residence". Explanatory variable "Area of residence" has 5 categories: Area1, Area2,...,Area5 and for the big school their coefficients are calculated to be $\hat{\beta_1}, \hat{\beta_2},...,\hat{\beta_5}$

I'm also modelling the students mood in a really small school, and we do not have much data. Some researcher says that for all the categorical explanatory variables we calculated before (for the big school), we can just use those calculated coefficients as restrictions in the analysis for the small school. E.g. in a lot of new statistical software one has the option to "save" a bunch of coefficients calculated for an explanatory variable, basically giving us a special function that can be used later in another GLM/GAM analysis (for the same categorical variables).

For the little school we a scarce amount of data such that that none of the categories (Area1-Area5) have significant p-values (we are a research group of "let's just evaluate p-values"). Using the restrictions we calculated for the big school, we have the model:

$log\frac{\pi}{1-\pi} = \beta_{new}[Area1=\hat{\beta_1}, Area2=\hat{\beta_2}, ... , Area5=\hat{\beta_5}]$

Only $\beta_{new}$ is estimated in the model for the little school, while the other betas are "restricted". The idea, or justification, is that the "Area of residence" variable affects the students mood exactly the same form (relationship between categories are the same) in both schools, except that the "effect" can be dampened or heightened depending on the MLE estimation of $\beta_{new}$.

Now, imagine if you do this with, say, 10 different variables and evaluate the p-values of the 10 coefficients $\beta_{new_1}, \beta_{new_2}, ..., \beta_{new_{10}}$. Some p-value in $\beta_{new_i}$ be significant due to chance and the wrong conclusion "the relationship between the categories of variable "x" is the same in the two schools" is drawn.

Question 1: Is this not another fancy version of data dredging from stepwise regression techniques (good answer here)?

Question 2: This is basically an attempt to be innovative, and be able to use information drawn from a big source, and extrapolate it towards a smaller source. Am I correct in thinking that this could be a good idea if one believed strong heartedly that e.g. "Area of residence" must behave the same in the two schools? But a disastrous idea when blindly fumbling in the dark trying to feel out p-values to determine which variables behave the same (have the same form) across the two schools?

Question 3:@Repmat points out in his answer that choosing the right functional is not critical. And, if I understand correctly, if it were critical, you would see it from your testing and validation sets. But what if the method described above would be used in the making of all the models (because it is common belief that it is a good method)? Then, would I not be comparing just bad models - leaving me with the least bad model?

Thoughts and reference request:Viewing this it made me think of how disastrous doing GLM/GAM analysis with the wrong functional of some covariate could be. For example, if one were able to fit $E[y] = x^2$, even though $E[y] = x$ were a more true model (like I tried to explain above), this would be horrible for future predictions of $E[y]$. Does any research exist on the consequences of choosing the wrong functional?

EDIT: fixed some of the structure of this question.

• Just on your last point - the range of your data is crucial here. If you are interpolating, then functional form doesn't matter so much - you'll get worse results (if you pick the wrong form), but still be "in the ballpark". Extrapolation is where is matters much more. Eg if you want to predict at $x=1000$ but your training data has values between $25$ and $100$ the choosing $x^2$ or $x$ will make a big difference. – probabilityislogic May 1 '16 at 13:33
• Thank you for pointing that out. Can I draw from your comment (which I think revolves around the form of a continous explanatory variable) that this goes for form of a categorical explanatory variable as well (basically what I am asking in the first two parts)? I feel like there is a big difference here... – Erosennin May 1 '16 at 13:48
• The equivalent for categorical variables is "new" categories. For example, if you had to predict for a school that had students in "area 6". – probabilityislogic May 1 '16 at 14:59

One way to think of what you're doing, that shows this is not data dredging. If you consider the big + little schools as one data set (person $i$, area $a$, school $s$), then you have a general model of:

$$\eta_{ias}=\beta_{as}$$

And your restricted model corresponds to the following assumption (b=big school, l=little school):

$$\frac {\beta_{ab}}{\beta_{al}}=\beta_{new}$$

And you can write the restricted model as :

$$\eta_{ias}=\beta_{ab}\beta_{new}^{I\{s=l\}}$$

You can thus fit into the framework of likelihood ratio testing, as your hypothesis consists of restricting the parameter space of the full model (the full model is where the little school gets each area effect estimated separately). This is one way to get an overall effect of whether the simplification is doing a good job.

However, this is a non-linear GLM, so the distribution under the null hypothesis may not be chisquare. Might be good to use parametric bootstrap or something (see below for rough outline)

This non-linear GLM can be iterated between two linear GLMs. Fix $\beta_{new}$ then estimate $\beta_{ab}$. Then fix $\beta_{ab}$ and estimate $\beta_{new}$. The method you've written is one cycle of this approach, setting the initial value for $\beta_{new}$ equal to zero. The justification for only using one iteration is that the big school has more data - so the influence of the little school data has on estimating the common parameters is negligible. The second iteration only involves the units in the small school.

Additionally, you don't need "restrictions" to estimate these models - just a modified design matrix for each regression. In the first model ($\beta_{new}$ fixed) instead of 0-1 dummy variables for each area, you have $\beta_{new}^{I\{s=l\}}$ instead of the "1" (i.e. replace the "1"s in the design matrix with $\beta_{new}$ for the small school units). For the second model, the design matrix is the linear predictor from the first model ($\beta_{ab}$) and the fitted model has no intercept.

But importantly, because the model is not strictly a GLM, you may not be able to rely on the standard normal approximation theory to get your p-values. You should check that this is the case. This is not hard to do:

1. simulate data from your model using the betas you estimated. Keep the lop sided big/little school samples, and keep the same area distributions (i.e. only change the binary response you are modelling)
2. Refit the model using the simulated data in step 1, getting simulated parameter estimates
3. Repeat steps 1&2 many times (say $500$ or so)

The simulated betas should have an approximate normal distribution. You should also check that the variance estimates for beta is close to the variance of the corresponding simulated betas. You can also refit the alternative, more flexible model, and check the likelihood ratio distribution in the simulations.

• Thank you very much! Wrt to data dredging, I was referring to the act of trying to fit many variables with restricted forms, not just one, without evaluating properly if the form should be the same in both schools. I'm not sure if that was clear in my question. Wrt estimation of $\beta_{new}$ I just thought this was done by regular MLE like any other case - that is - when separating the big school and the little school. Anyway, I appreciate the explanation of the estimation procedure! Thank you! – Erosennin May 1 '16 at 15:40
• should the definition of $\beta_{new}$ be the inverse of the fraction above? also, what you mean by non-linear GLM confuses me a bit to be honest – Erosennin May 1 '16 at 16:30
• You are right - it should be inverted. With regards to non-linear GLM, you can't write the predictor in the form $\eta_i=x_i^T\beta$ for some known vector $x_i$ and parameter vector $\beta$. In terms of the hypothesis testing - the likelihood ratio test can incorporate testing multiple constraints simultaneously. So you are not doing multiple tests - just one. If it fails, then you can look at which components are driving the failed test. Not data dredging, but providing some insight into a large and complex test (by breaking into smaller pieces). – probabilityislogic May 1 '16 at 22:25

You are asking the central question for any parametric type of analysis, what if my model (in casu the functional form) is wrong? This is not as depressing as one might think, you lose unbiasedness but can maintain consistency provided the regressors are truly exogenous. Therefore the relevant question is really if you have omitted variables. If you do the estimates are simply wrong, and nothing (in general) can bail you of that assuming you don't have an instrument.

You mention predictions, if this is your interest then consistency might be less important and the above should not worry you to much. Your test and validation sets will tell you, if you are doing well enough.

In any case there is little you can do, other than grab more data or make do with what you got.

• Thank you for answering! I was actually just about to edit the question, but hope you understood what I tried to ask without too much effort! I'll try to follow up with a couple of questions here, if you don't mind: 1) wrt predictions, of course you are right, a test and validation set will clear things up. What bothers me is the method used. It feels to me like what @Frank Harell in his book talks against: data dredging and trusting p-values. 2) wrt wrong functional form not being critical, you are probably correct. ..ctd.. – Erosennin Apr 30 '16 at 20:20
• ...ctd.. but I'm just worried that by trusting a form by it's p-value (of $\beta_{new}$) is dangerous. What if the first half of the coefficients, say $\beta_1,\beta_2, \beta_3$ fit well to the new data (little school in my example), but $\beta_4$ and $\beta_5$ fits just horribly! The p-value of $\beta_{new}$ might be really small. I'm trying to make do with what I got wrt little data, I did not quite understand what you wrote about exogenous variables, though, but found a nice explanation (stats.stackexchange.com/questions/59588/…) – Erosennin Apr 30 '16 at 20:29
• You can never test if a given form is right or not, but what you do is measure if one form performs better than another given some relevant criteria. Ofc, you should never blindly trust a pvalue, instead try and determine what the model should be able to do and see how well it works. – Repmat Apr 30 '16 at 20:34