Interpretation of Linear Regression Models That Include Transformations or Interaction Terms W. Dana Flanders, MD, ScD, and David S. Freedman, PhD
Rebecca
DerSimonian,
PhD,
In linear regression analyses, we must often transform the dependent variable to meet the statistical assumptions of normality, variance stability, or linearity. Transformations, however, can complicate the interpretation of results because they change the scak on which the dependent variable is measured. In this setting, the inclusion of product terms or the transformation of some independent (or predictor) variables may further complicate inrerpretarion. In this artick, we present some interpretations of linear models that include transformations or product terms. We illustrate these interpretations using regression analyses designed to study determinanrs of serum testosterone levels. These examples show how one can present results using simpk measures, such as medians, and interpret regression parameters. Ann Epidemiol I992;2:735-744. Epidemiology, epidemiologicalmethods, statistical models, regression.
KEY WORDS:
INTRODUCTION In ordinary linear regression analyses, we must often transform the dependent to satisfy the usual statistical assumptions of normality,
variable
variance stability, and linearity
(1, 2). For example, many biologic variables have an approximate log-normal distribution so that the logarithmic
transformation
may be appropriate. The use of transforma-
tions, however, complicates the interpretation of model parameters. Such difficulties arise because the transformation changes the scale on which the dependent variable is measured. The inclusion of product terms or the transformation of some independent variables may further complicate
interpretation
in these models.
In this article, we present some simple interpretations transformations transformation original,
of linear models that include
or product terms. In particular, we show how use of the inverse leads to a model for the median of the dependent variable on the
untransformed
more traditional to understanding
scale. We do not suggest that these interpretations
presentations,
but that they complement
replace
and enhance them as an aid
results.
NOTATION Let Y denote
the dependent
variable of interest;
2 = T(Y)
denote
the dependent
variable after transformation by T; and X, , X,, . . . , X, denote the independent (predictor) variables. We assume the following linear model:
E(Z) = EUO’)) = P(Z) = & +
p,g(x,)
+ . . . +
pexe,
(1)
From Emory University School of Public Health, Atlanta, (W.D.F.); National Institute of Child Health and Human Development, Bethesda, MD (R.D.); and Centers for Disease Control, Atlanta, (D.S.F.), GA. Address reprint requests to: W. Dana Flanders, MD, ScD, Emory University School of Public Health, Division of Epidemiology, 1599 Clifton Road, Atlanta, GA 30329. Received March 27, 1991; revised June 21, 1991. 0 1992 Elsevw Science Publishing Co., Inc.
1047~2797/92/$05.00
Flanders et al. INTERPRETATION OF REGRESSION MODELS
TABLE 1
Multiple
AEP Vol. 2, No. 5 September 1992: 735-744
linear regression results using square root of serum testosterone
levels
Predictor Intercept Log BMI (kglm2)
68.39
1.60
-9.71 -0.29 0.53
Age (Y) Race (white = 0, black = 1)
0.41 0.03 0.19
SE = standard emx estimate; BMI = body mass index
where E(Z) denotes the expected value transformation, g, of X, to emphasize that
or mean of 2. We explicitly include the independent as well as dependent variables
may require transformation (1). Moreover, a predictor variable, denote an interaction or a product term, such as X2 * X,.
say X,,
could
also
models
that
Example To illustrate include
the difficulties
transformations
that
or product
may arise with terms,
consider
interpretation the analyses
of linear conducted
and colleagues (3) of data from the Centers for Disease Control veterans (4-6). They studied determinants of serum testosterone linear regression models, including:
by Freedman
study of US Army levels using several
E(Z) = Po + P, lodx,) + &XL + PJ,, where
2 is the square
root of testosterone
levels
(ng/mL);
(2)
X, is the body mass index
(BMI; kg/m’); X2 is age (years) at examination; and X, is race (0 for whites and 1 for nonwhites). As summarized in Table 1, the parameter estimate for log(X,) is -9.71, which,
by the usual interpretation,
indicates
that
the square
root of the testosterone
level decreases, on average, by 9.71 per unit increase in the logarithm of the BMI. The practical meaning of this result is unclear, in part because the transformations alter or distort
the usual, arithmetic
scale for measurement
of testosterone
BMI. Researchers and practitioners may be unaccustomed level or BMI on the square root or logarithmic scale.
INTERPRETING
MODEL-PREDICTED
concentration
to thinking
and
of testosterone
MEANS
Let the transformation T(Y) of equation 1 be in the family of power and denote the inverse transformation by T-‘(Z). Thus:
2 = T(Y) = Y”
and T-‘(Z)
= 2 “*
transformations,
ifh # 0, (3)
= log(Y)
and T-‘(Z)
= exp(2)
ifh = 0.
The definition of T(Y) as the natural logarithm for h = 0 is a common notational convention that we adopt here. When the transformed variable 2 = T(Y) has a gaussian distribution, the results in Appendix 1 show that we can estimate the median of Y by:
AEP Vol. 2, No. 5 September1992: 735-744
Flanderset al. INTERPRETATION OF REGRESSION MODELS
Q = T-l@)
= fill”
=
= exp(p)
tp,
+
= exp(&
737
&g(X,) + * . . +&x,)“”
ifh # 0,
+ &g(X,)
ifX = 0,
+ . . . + &X,)
(4)
where fi is the estimate of the mean of 2.
INTERPRETATION
OF M AS A GENERALIZED
When A = 0 (the logarithmic of Y, so one can interpret
MEAN
transformation),
T-‘(P)
fi estimates the mean of the logarithm
as an estimate
of the geometric
mean of Y (7)
(Appendix 2). In this context, we can regard 0 = T-‘(P) as a generalized mean estimate of Y where fi is the model predicted mean of T(Y). In particular, we can interpret M as an estimate of: the geometric h = - 1, or the arithmetic mean if A = 1.
Example
mean if X = 0, the harmonic
mean if
(continued)
To illustrate the use of equation 4, consider further the data analyzed by Freedman and colleagues using the model given by equation 2. Equation 4 and the results in Table
limply that the median testosterone
level is M = (68.39
-
9.71 * log(26)
-
0.29 * 37 + 0.53 * O}’ = 677, for a white veteran of average age (37 years) and BMI (26 kg/m2); and is &I = (68.39 - 9.71 * log(26) - 0.29 * 37 + 0.53 * 1}2 = 705, for a similar black veteran. Alternatively, as indicated above we can also interpret these values as estimates of a generalized mean for white and black veterans, respectively.
INTERPRETING INDEPENDENT
THE ASSOCIATION VARIABLES
OF THE
DEPENDENT
WITH
THE
When, in addition to the transformation of the dependent variable, the model includes a transformation of an important independent variable, further difficulties arise in interpreting
regression results. Suppose that interest centers on the association between
Y and X, in the model given by equation 1. One usually interprets p, , the parameter associated with X1, as implying that the mean of T(Y) changes by p, per unit increase in g(X,). This “usual” interpretation is precise, but the association of Y and X, is not intuitively clear since it involves two transformations, T(Y) and g(X,). However, we can estimate the effect of Xi on the median of Y by taking derivatives on both sides of equation 4:
6!ibx, = (8, + &g(X,) = exp(P,
+ . . . + P~X,)(l’“-L’~,g’(X,)/A
+ &(X1)
+ . . . + &XJPg’W,)
ifh # 0, ifh = 0,
(5)
where g’(X,) is the first derivative of g with respect to X, , assumed to be an interval n variable. This expression for 6M/6X, estimates the slope: the change in the median of Y per unit change in X,,. Use of 6&I/6X, provides us with an alternative way to understand the effect of X1 on Y. The slope characterizes the association between Y and Xi, yet avoids the possibly unfamiliar, transformed scales of measurement. We note that the association can depend on the covariates, so one may need to illustrate the association between Y and
738
Flanders et al. INTERPRETATION
X,,
OF REGRESSION
for several combinations
AEP Vol. 2, No. 5 September 1992: 735-744
MODELS
of values. Furthermore,
when the model includes no
transformations, equation 5 reduces to the usual interpretation of the coefficients. This follows mathematically since no transformation implies X = 1, g’(X) = 1 and equation 5 reduces to 6&l/6X,
= B,.
Example (continued) To illustrate the use of equation
5, we again consider the results in Table
1. For a
black veteran (37 years old, BMI of 26) the model predicts that the median testosterone level would change by: &/6X,
= -2
* [68.39 - 9.71 * log(26)
- 0.29 * 37 + 0.53 * 1][9.71 * l/26] (6)
= - 19.8 ng/dL per unit change in BMI,
where we have usedg(X1) = log(BMI), andg’(X,) = l/BMI. Thus, the model predicts that for black males, the median testosterone level would decrease by about 20 ng/dL per unit change in the BMI,
INTERPRETATION
at the average values of BMI and other predictors.
IN THE PRESENCE OF EFFECT MODIFICATION
The preceding
results yield an expression
for &I as a function
of the independent
variables; the expression is stratum-specific in that it pertains to a specific combination of the X’s. We can use these results to express the generalized mean of Y for a standard population. Use of a standardized mean may be particularly helpful if interaction or product terms have been included in multiple linear regression models, since such terms may further complicate the presentation of summary results. To illustrate, consider the model of equation 1 where X, = X1 * X, . Here we focus on the association between Y and X, . This model implies that the change in the mean of T(Y) per unit increase in X, is (& + &,X3). This, in turn, implies that the effect of X, depends on or is modified by X, . Although many situations may require detailed presentation of the complete pattern of effect modification, mean, may complement
a summary measure, such as a standardized
that presentation.
Standardization is a frequently used epidemiologic technique that permits the epidemiologist to study effects of the variable of interest, say X, , after controlling for confounding (8-12). Standardization yields a summary measure that is a weighted average of effects in which the weights reflect the distribution of the covariates in a standard population.
The standardized measure is interpretable
as an average effect
and can complement the presentation of complex patterns of effects. Greenland (9) argued that standardized estimates are meaningful summary measures since they are interpretable as average values for the standard population. For instance, if Xlk, X,, , X,,,. . .k=l,... N are the covariate values for individual k in the standard population and g(X,), x,, x,, . . . are the corresponding standard population, then for the model at hand:
estimates the mean of T(Y)
for the standard population,
average values in the
as a function
of g(X,),
X, ,
AEP Vol. 2, No. 5 September 1992: 735-744
TABLE 2 Multiple linear regression results using square root of serum testosterone Model that includes a product term Predictor
739
Flanders et al. OF REGRESSION MODELS
INTERPRETATION
levels: SE
B
1.664 0.41 0.03 2.83 0.07
67.34 -9.70 -0.26 8.28 -0.21
Intercept Log BMI Age Race Race * age“ a race times age “product” term. SE = standard error estimate; BMI = body mass index.
x, . . . . Analogous population k,
to equation
4, the generalized
mean of Y for this standard
is:
= T-‘(&)
..-
= (& + &g(XJ
+ ,&zx, + . . . ifX # 0,
+ pp-l%~-* + up&)“” = exp(&
I+ prg(X,)
+ B,X,
+ . . .
+ Pp-&-1
ifX = 0.
f &px,x,)
(8) Moreover,
one can estimate the association
and X, for the standard population 6ti,/6X,
A-
= (p, + &g(X,)
between the generalized mean of Y
by:
+ . . . + &,X0_,
+ pp?qxy-qpz = exp(&
n+ p,g(X,)
+ @Q/h
ifA # 0, (9)
+ . . . + &_,Xp_,
+ &x,x,)
csz + P&J
ifh = 0.
In this expression, 6&l/8X, estimates the change in the (generalized) mean of Y for the standard population per unit change in Xz and we assume that the change in X, occurs uniformly remain fixed.
in the standard population
Example (continued) To illustrate the use of equations
and that other independent
variables
8 and 9, consider additional analyses by Freedman
and colleagues (3) which suggested that the association between testosterone level and age depended on race. Table 2 summarizes the parameter estimates from the following model:
E(Z) = P, + P,dX,) + PJz + PA + P&,
(10)
where 2, X, , X1, and X, are defined in equation 2, and X4 is the product or interaction term of age with race. To illustrate the implications of this model, presentation of standardized results is helpful. Thus, consider a population, like the Centers for Disease Control study population, in which the average age is 37.8 years, 12% are black, the average of the logarithm of BMI is 3.28, and th e average of the age * race term is 4.64. For this standard population, equation 8 shows that, assuming normality, the generalized mean testosterone level is:
740
Flanders et al. INTERPRETATION
TABLE
3
OF REGRESSION
AEP Vol. 2, No. 5 September 1992: 735-744
MODELS
Multiple linear regression results using logarithm of serum testosterone
Predictor
9.832 -0.768 -0.021 0.039
Intercept Log BMI Age Race
% Change in median testosterone level/unit increase in predictor variable
$E
B
levels
0.13 0.03 0.002 0.02
-53.6% - 2.1% 4.0%
SE = standard error estimates; BMI = body mass index.
~;i,=[67.34-9.7O*3.28-O.26*37.8+8.28*O.12-O.21*4.64]2
(11)
= 664. For this population,
equation
9 yields:
6~l6X2=2~[67.34-9.7O*3.28-O.26*37.8+8.28*O.12-O.21*4.64 (I12)
. [- 0.26 - 0.21 * 0.121 = - 14.6. Thus the model predicts that the median testosterone 15 ng/dL per year as the standard population
level would decrease by about
ages, at the standard values of age and
other predictors.
SPECIAL CASE: LOGARITHMIC DEPENDENT VARIABLE
TRANSFORMATION
OF THE
When the logarithmic transformation is used to transform the dependent variable (A = 0), some additional, simple interpretations are possible. Here, we also assume that product terms are not required (& = 0). For the logarithmic transformation, the usual interpretation
of /!?,is that it repre-
sents the change in the logarithm of Y per unit increase in g(X,). Using equation the estimated proportional change in the median of Y per unit change in X, is:
@mx,)/&l = &g’(XJ. Use of derivatives,
as in equation
5,
(13)
13, is inappropriate with a categorical
independent
variable. To derive an estimate of the proportional difference in the generalized mean for a categorical variable, we use equation 14 along the lines suggested by Stryker and coauthors (13). For example, to compare the mean among those with X2 = x1 with that among those with X, = x, (13), we have:
[ticc,=x,j - ~~~~=~,~l~~ji(~,=~,,~ = x,,)) exp&(x,
Example
-
1.
(14)
(continued)
To illustrate use of equation 13, consider again the data analyzed by Freedman and colleagues (3). Table 3 summarizes the parameter estimates for the model given by equation 2, except that 2 is the logarithm of testosterone level instead of the square root. The p for age is -0.021, indicating that the mean of the logarithm of testosterone
Flanderset al. INTERPRETATION OF REGRESSION MODELS
AEP Vol. 2, No. 5 September 1992: 735-744
741
TABLE 4 Multiple linear regression results using logarithm of serum testosterone levels: Model that includes a product term % Change in median testosterone level/unit Predictor
s
Intercept Log BMI Age Race Race * age
9.755 -0.768 -0.021 0.610 -0.015
E
increase in predictor” variable
0.132 0.033 0.002 0.227 0.006
53.6% -2.3% 4.4%
’All percentage changes are for a standard population in which the average of log(BM1) is 3.28, 12% are black, the average age is 37.8 y, and the average of the race by age variable is 4.64. SE = standard emx estimate; BMI = body mass index.
level decreases
by 0.021
given by equation is 100% (-0.021)
per l-year
increase
in age. An alternative
interpretation,
13, is that the proportional change in the median testosterone level = -2.1% per l-year increase in age. The right column in Table 3
indicates the percentage change in the median testosterone each independent variable.
level per unit change in
For logarithmic transformations of the dependent variable, the proportional change in the generalized mean of Y can be estimated even if the model implies effect modification.
In this instance,
the estimated (generalized)
mean for Y in the standard
population is exp(@,), where fi, is given by equation 7. Now the proportional change n in M for each unit increase in X, (for all members of the standard population) is, from equations
8 and 9,
where &ls is given by equation 8 with covariate values Xlk, X,, + 1, Xjk . . . for k = 1, 2, . . . N. Freedman
and associates
(3) also modeled the logarithm
of testosterone
so that
the effect of age depended on race. Table 4 summarizes the parameter estimates from the following model: E(Z)
=
P,
+
PJ,
where 2 is the logarithm of testosterone,
+
&X2
+
P3X3
+
and X, , . . . , X4 are the logarithm of BMI,
race, age, and the product term of race with age as in equation implications
(16)
/34x49
of this model, the presentation
10. To illustrate the
of standardized results is helpful. Thus,
consider a population, like the Centers for Disease Control study population (4-6), in which the average age is 37.8, 12% are black, the average of the race by age variable is 4.64, and the average of log(BMI) is 3.28. For a population like this, equation 15 predicts that the proportional the generalized mean testosterone level is 100 * (-0.021 - 0.015 * 0.12) per year due to aging of the standard population.
change in = -2.3%
DISCUSSION We have illustrated, by example, how interpretation may be difficult if one uses transformations in linear regression. Another illustration is provided by Roidt and coworkers (14), who used responses on a food frequency questionnaire to predict serum
742
AEP Vol. 2, No. 5 September 1992: 735-744
Flanders et al. INTERPRETATION OF REGRESSION MODELS
beta-carotene
and
retinol
levels.
To simplify
interpretation,
transformations to calculate or present means distributions were skewed, so they transformed
they
did not
use the
and medians. On the other hand, variables before calculating P values
and before doing linear regression analyses. In this work, Roidt and coworkers faced the difficulties associated with the use of transformations, and simplified the presentation by not
using
technically
transformations acceptable,
for calculating
means
and medians.
approach
is
were transformed
in
some analyses but not in others, and that regression parameters terms of statistical significance but not in terms of the meaning
were interpreted and implications
in on
untransformed
scale.
In this
that
Their
variables
a familiar,
but has the disadvantages
article,
we discussed
and
illustrated
some
interpretations of linear models that may simplify presentation when the models include transformed variables like those encountered by Roidt and associates. The key interpretations
follow
because
use of the inverse
transformation
allows one to model
the median of the dependent equation 4. The investigator
variable on the original, untransformed scale as given by can calculate approximate standard errors using the delta
method (15). Some of the approaches
presented
example,
Stryker
between
beta-crotene
regression
models,
and
colleagues and smoking,
they
here were used in the published
(13)
presented
their
estimates
with use of the format
predicted
the mean
For
association
here.
From linear
described
of the logarithm
literature.
of the
of beta-carotene
levels
using the number of cigarettes smoked per day and other variables. The estimated cient for cigarettes (packs/d) was -0.33. Thus, they estimated that the geometric
coeffimean
beta-carotene levels for men who smoked 1 pack per day was exp( -0.33) = 72% of the levels for nonsmokers. Their interpretation resembles that presented here: Equation 14 yields a proportional consistent
difference
with the interpretation
of lOO%[exp(-0.33) used by Stryker
formation of the dependent variable metric mean estimate are equivalent.
-11
= -28%.
and associates,
(h = O), the model predicted Centers for Disease Control
however, mean,
approximate because
the
arithmetic
properties
are well documented
mean.
This
of some means, (7). Another
limitation reflects
(4-6). mean
is mitigated
such as the harmonic
limitation
is
median and the georesearchers also used
similar interpretations in presenting results of their studies of veterans The approach has some limitations. For example, a generalized closely
This result
since for the log trans-
need
not
somewhat,
or the geometric
the potential
dependence
of the slope on covariate values or their distribution in the standard population. Although some of the difficulties that arise from the use of transformations might be avoided by using nonlinear regression or generalized linear models (16), linear regression
analyses
and
transformations
are commonly
used so that
the present
approach
remains relevant. Although most of the ideas presented here are not new, they are not widely used in the epidemiologic literature. They provide a simple way to explain results of regression analyses and could supplement usual modes of presentation. A key advantage is that use of these ideas allows the epidemiologist to present results on the original, untransformed scale that is familiar to clinicians, rather than on a transformed scale that involves logarithms or power transformations.
APPENDIX
1 with mean Assume that the transformed variable 2 = T(Y) h as a Gaussian distribution p and variance oz. Since such a distribution is symmetric about the mean,
AEP Vol. 2, No. 5 September 1992: 735-744
INTERPRETATION
0.5 = Making gives:
the substitution
/:m
(1/27rf_&}i/2exp(-(7
T(Y)
- ~)~/2t9)dz.
(la)
= 2 and applying the change-of-variable
0.5 = 1,‘-11”, {1/27&}“2
exp(-(T(y)
743
Flanders et al. OF REGRESSION MODELS
technique
- ~)~/2a~)T’(y)dy, (2a)
=I ,T_l;:;,f(r)dr> where T’(Y) is the first derivative of T with respect to Y and f(Y) is the density function of Y. Equation la follows from the symmetry of the Gaussian distribution of T(Y)
around its mean p. Equation
T-‘(p),
so by definition T-‘(p)
2a shows that half the values of Y are less than
= M is the median of Y. In other words, the median
of the distribution
of Y equals the inverse transformation of /L, the mean of T(Y), where @ is provided T(Y) is G aussian; we can estimate the median of Y by T-‘(p), the estimate of the mean of T(Y). In practice, T(Y) may only be approximately Gaussian, making the estimates approximate. For the power transformations, these results show that one can estimate the median of Y as:
fi =
T-~(P)
=
p”A
(& + fi,g(X,) + . . . + gpxpy*
=
= exp(/L) = exp(j& + b,g(X,)
APPENDIX
+ . . . + &X,)
ifh # 0, ifh = 0,
(3a)
2 As suggested elsewhere
(7),
the arithmetic,
harmonic,
and geometric
means
are
included in the family: M = T-‘(C G(Y)
T(Y,)lN)
= T-‘E(T(Y))
for the sample mean, or
(4a)
for the population mean,
(5a)
where the summation in equation 4a is over the Y, in the sample and T(Y) is given by equation 3. We substitute the estimated value of E(T(Y)) into equation 5a to estimate the generalized mean. To obtain the generalized mean for a standard population, for each combination
of independent
we average the mean
variables in the standard population
and then
retransform:
M, =
T-l
(64
where the summation is over the N individuals in the standard population and E(T(Y,)) denotes the expected value of T(Y) for covariate values Xlk, X21,, XJk, . . . . Equation 8 follows by substituting model-predicted estimates of the means into equation 6a.
744
Flanders et al. INTERPRETATION
OF REGRESSION
AEP Vol. 2, No. 5 September 1992: 735-744
MODELS
REFERENCES 1.
Kleinbaum
Multivariable
DG,
Methods.
Kupper
LL, Muller
2nd ed. Boston:
KE. Applied
PWS-Kent;
Regression
Analysis
and Other
1988.
2.
Box GEP,
Cox DR. An analysis of transformations,
3.
Freedman
DS, O’Brien
J R Stat Sot Ser (B).
1964;26:
211-43. of testosterone Arterioscl
Thromb
4.
Centers
Psychosocial 5. Physical
Centers
7.
for Disease Control
JJ. Relationship in men,
Vietnam
experience
study. Health
status of veterans.
1:
Vietnam
experience
study. Health status of veterans.
experience
study. Health status of veterans.
II:
JAMA.
Vietnam
III:
1988;259:2715-9.
Geometric
mean.
In: Encyclopedia
of Statistical
Sciences.
v. 3. New
1983.
8.
Freeman
9.
Greenland
Stat Med.
F, Barboriak
and other characteristics
1988;259:2701-7.
for Disease Control
outcomes, Wiley;
DeStefano
cholesterol
1988;259:2708-14.
Hines WGS.
York: John
JAMA.
for Disease Control
JAMA.
Centers
Reproductive
Flanders WD, lipoprotein
1991;11:307-15.
characteristics,
health,
6.
TR,
levels to high-density
DH, Holford TR.
Summary
S. Interpretation
rates, Biometrics.
and estimation
1980;36:195-205.
of summary ratios under Heterogeneity,
1982;1:217-27.
10.
Lane PW,
prediction.
Nelder
Biometrics.
11.
Miettinen
OS.
12.
Wilcosky
TC,
adjustment 13.
Stryker
Analysis
Standardization Chambless
of epidemiologic
and alcohol
JA.
of covariance
and standardization
as mstances
of
1982;38:613-21.
WS,
of risk ratios,
LE. A comparison
measures,
J Chronic
Kaplan LA, Stein
consumption
Am J Epidemiol. of direct
and regression
Dis. 1985;38:849-56.
EA, et al. The relation
to plasma beta-carotene
1972;96:383-8.
adjustment
of diet, cigarette
and alpha-tocopherol
smoking,
levels, Am J Epidemiol.
1988;127:283-96. 14. estimates 15. MA: MIT 16. and Hall;
Roidt L, White
E, Goodman
GE, et al. Association
of food frequency
of vitamin A intake with serum vitamin A levels, Am J Epidemiol. Bishop YMM,
Fienberg SE, Holland PW. D’lscrete Multivariate
questionnaire
1988;128:645-54.
Analysis.
Cambridge,
Press; 1975. McCullagh 1989.
I’, Nelder JA. Generalized
Linear Models.
2nd ed. New York: Chapman