c / i/ PREVENTIVE VETERINARY ...
c
/
i/
PREVENTIVE
VETERINARY
MFD1cINE
ET ‘?F;VTI;R
Pwventive Velerinarv
I_IU”U
. .&.a*.
_
.-.-
_.__.-
Medicine 46
__......-,
(2000) 225-247www.ekevier.nl/locate/prevetmed
Graphical approaches to support the analysis of
linear-multilevel models of latib pre-weaning
growth in Kolda (Senegal)
Renaud Lancelota7*, Matthieu LesnoffQ, Emmanuel Tillardb,
John J. McDerrnottc
%Zentre de Coopération Internationale en Recherche Agmnomique pour le Développement, Département
Élevage et Médecine Vétérinaire (CIRAD-EMVT),
TA 3O/A, Campus International de Baillarguet, 34398
Montpellier Cedex 5, France
bCIRAD-Eh4Vr Station de Ligne-Paradis, 7 chemin de I’IRAT, Ligne-Paradis, 97410
St Pierre, Ile de la Réunion, France
cIntemationai Livestock Research Institute (ILRI), PO Box 30709, Nairobi, Kenya
Received 28 December 1999; accepted 30 May 2000
Abstract
Linear-multilevel models (LMM) are mixed-effects models in which several levels of grouping
may be specified (village, herd, animal, . . .). This study highlighted the usefulness of graphical
methocis in their analysis through: (1) the choice of the fïxed and random effects and their structure,
(2) the assessment of goodness-of-fit and (3) distributional assumptions for random effects and
residuals.
An LMM was developed to study the effect of ewe deworming with morantel on lamb pre-
weaning growth in a field experiment involving 182 lambs in 45 herds and 10 villages in Kolda,
Senegal. Growth was described as a quadratic polynomial of age. Other covariates were sex, litter-
size and treatment. The choice of fïxed and random effects relied on three graphs: (1) a trellis
display of mean live-weight vs. age, to Select main effects and interactions (fixed effects); (2) a
trellis display of individual growth curves, to decide which growth-curve terms should be included
as random effects and (3) a scatter plot of parameters of lamb-specifk regressions (live-weight vs.
quadratic polynomial of age) to choose the random-effects covariance structure.
Age, litter-size, agexlitter-size,
litter-size x treatment and age xlitter-size xtreatment were
selected graphically as fïxed effects and were signifïcant (~~0.05) in subsequent statistical
models. The selection of random-effect structures was guided by graphkal assessment and
* Corresponding author. Present address: Institut Sénégalais de Recherches Agricoles, Laboratoire National
d’Élevage et de Recherche Vétérinaire (ISRA-LNERV), BP 2057, Dakar-Hann, Sénégal. Tel.: +221-832-49-02;
fax: +221-821-18-79.
E-muil uddress: renaud.lancelot@cirad.fr (Ii. Lancelot).
0167-5877/00/$ - see front matter 0 2000 Elsevier Science B.V. AL1 rights reserved.
PII: SO167-5877(00)00155-0
226
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
comparison of the Akaike’s information criterion for different models. The final random-effects
selected included no random effect at tbe village level but intercept, age and squared-age at the
herd and lamb levels. Tbe structure of tbe random-effects variance-covariance matrices were
blocked-diagonal at the herd level and unstructured at the lamb level. An order-1 autoregressive
structure was retained to account for serial correlations of residuals. Smaller residual variante at
90 days than at younger ages was modeled witb a dummy variable taking a value of 1 at 90 days
and 0 eIsewhere.
Ewe-deworming with morantel during tbe rainy season lead to higber latnb live-weights
(probably related to a better ewe-nutrition and -health status). A positive correlation was
demonstrated between early weight and growth rate at the population level (with important lamb
and herd-level random deviations). Tbe persistence of .this correlation at older ages should be
checked to determine whether early weigbts are good predictors of mature weights and ewe-
reproductive lifetime performance. 0 2000 Elsevier Science B.V. Al1 rights reserved.
Keyworuk: Lmear-multievel growth model; Random effect; Grdphical methods; hz-weaning growth, Castro-
intestins parasitism; Mordntel; Sheep; Senegd
1. Introduction
Studies in veterinary medicine frequently measure phenomena repeated over time
(such as live-weights of the same animais). The analysis of these measures as multivariate
outcomes over time is usually more informative than repeated cross-sectjonal studies, e.g.
live-weights at given ages or cumulated production data such as 305day milk yield for
dairy cows (Grohn et al., 1999). In addition, veterinary studies ofien rely on multilevel
designs, in which a number of units at different levels (for example: villages, herds within
villages, and animals within herds) are sampled. Data arising from these studies often
have both longitudinal and multilevel features. Mixed-effect models have been used to
analyze such longitudinal and multilevel data appropriately (Diggle et al., 1994). They
provide a choice of effect type (either fixed or random) depending on study design and
objectives at different grouping levels. The statistical properties of linear-multilevel
models have been well described in recent years (e.g. Goldstein, 1986, 1995’; Pinheiro and
Bates, 1996) and different software developed to fit them (Kreft et al., 1994; Littel et al.,
1996).
Linear-multilevel growth models have been applied recently in veterinary studies (e.g.
Puyalto et al., 1997; Green et al., 1998). However, because of their complexity and the
different inferences that cari result from different constructs of fixed and random effects
(McDermott et al., 1997), tare must be taken in mode1 specification. Model-specification
issues include the choice of fixed and random effects at different levels and the structure
of random-effect variance-covariance matrices. The recent availabjlity of software
combining computing tools and graphjcal methods, provides new and intuitive tools for
better exploring the most-appropriate mode1 structures prior to fitting and in assessing
and comparing fitted models (Pinheiro and Bates, 2000).
In thjs paper, we emphasize the use of graphical methods as a tool in mode1 building
and assessment for repeated and multilevel outcomes. As an example, a linear-multilevel
R. Lancelot et al. /Preventive Veterinmy Medicine 46 (2000) 225-247
227
growth mode1 was developed to study the effect of ewe deworming with morantel on
lamb pre-weaning growth in a field experiment in several herds and villages in Kolda,
Senegal. Graphical and statistical analyses focused on the following objectives: (1) to
estimate the effect of morantel treatment on lamb growth correctly, (2) to derive
population-averaged estimates for lamb growth curves, (3) to evaluate the importance of
lamb, herd and village-level deviations around the populations values and (4) to provide
information on factors influencing live-weight growth in individual lambs under field
conditions. The advantages of including graphical methods in mode1 building and
assessment are highlighted.
2. Materials and methods
2.1. Data
Kolda is located in the sub-humid, southem zone of Senegal. From 1983 to 1995, mean
annual rainfall was 963 mm (mainly occurring between June and November). Means of
maximum- and minimum-daily temperatures were 35.5 and 20.6”C, respectively. The
landscape was a mosaic of forested plateaus and lowlands cultivated with rice, cotton,
maize, and groundnuts. A detailed description of the farming systems in this area was
given by Faugère et al. (1990). Briefly, the main ethnie group was the Fulani. Crop
farming was the dominant activity but livestock production was widely practiced (cattle,
sheep and goats). Sheep were mainly Djallonke (a trypanotolerant West African dwarf
breed reared extensively for meat). Normally, they were left to graze freely in the bush
except during the rainy season when they were tied to avoid trop damage. A high
incidence of gastrointestinal helminthoses occurred in sheep during the rainy season
(Tillard, 1991). An abattoir survey revealed that Trichostrongylus colubrifomis,
Oesophagostomum columbianum, Haemonchus contortus and Strongyloides papillosus
were the most-prevalent parasites in this region (Fritsche et al., 1993).
A field experiment was designed to assess the profitability of deworming sheep flocks
in this zone and was conducted from July 1987 to June 1988. Partial results were
published previously (Faugère et al., 1990; Tillard, 1991). Ten villages were selected
according to agro-ecosystem representativeness and accessibility criteria. Herds
belonging to volunteer farmers were enrolled and all their animals ear-tagged (Faugère
and Faugère, 1986). Trained surveyors recorded demographic events and growth data
during fortnightly visits. Al1 10 villages enrolled in the study were randomly allocated to
one of two treatment groups (morantel or control) using a random-number table. Five
villages were allocated to the morantel group, and five villages to the control group. Forty-
five herds were included in the experiment: 17 in the control villages and 28 in the
morantel villages. Treatments were applied at the village level. In the morantel-treated
villages, ail sheep older than 3 months were drenched three times during the rainy season
with an oral administration of morantel (Exhelm@, Pfizer) at 7.5 mg kg-’ live-weight.
No sheep in the control villages and no lamb younger than 3 months in the morantel
villages were given any dewormer. Thus, the treatment-effect on lambs was mediated
through its effect on their ewe.
228
R. Lancelot et a[. /Preventive Veterinaty Medicine 46 (2000) 225-247
In Kolda, lambing distribution is bimodal, with two peaks from October to December
(end of the rainy season) and March to May (end of the hot, dry season) (Clément et al.,
1997). Tbese two lambing seasons define two cohorts of lambs subject to different
nutritional and environmental conditions under different husbandry practices (Faugère
et al., 1990), all with possible important effects on live-weight growth. Because the
purpose of this paper was mainly methodological, we decided to limit the study to the
cohort of lambs born from 1 October to 31 December 1987. A more-comprehensive study
of the whole data set Will be presented in subsequent papers.
There were 182 lambs with 922 weight measurements in this data subset. At the time of
data collection (1987-1988), follow-up information was checked and stored in temporary
files (Faugère and Faugère, 1993). Weights stored in these files had been adjusted for the
standard ages of 15, 30,45, 60,715 and 90 days (raw weights were thrown away). Before
analysis, these files were gathered in a relational database (Lancelot et al., 1998). Given
that each lamb was weighed every 2 weeks, little information was lost on the curvature of
live-weight-growth trajectories. However, because of this standardization, a few weights
(6%) were discarded. We believe this loss of information had no important effect on our
results, beyond a small decrease in statistical power.
Previous studies suggested that sex and litter-size (number of lambs bom) were
important factors to consider for lamb pre-weaning growth (Poivey et al., 1982; Fall et al.,
1983; Tuah and Baah, 1985; Armbruster et al., 1991; Abassa et al., 1992). T~US, sex
(male or female), litter-size (single or multiple: twins and triplets) and treatment
(morantel or control) were the flxed-effect classification variables used.
2.2. Data analysis
2.2-l. Craphical analyses
Graphical analyses followed the three-step exploratory data-analysis (EDA) framework
proposed by Cleveland (1993): (1) display the data with respect to its inherent structure
and the mode1 to be fitted; (2) display the fitted mode1 and check its assumptions and (3)
examine mode1 consistency with the original data. TO enhance visual comparisons of
different variable combinations within hierarchically structured data, Cleveland (1994)
proposed the method of trellis graphies for which software was subsequently developed
(Becker et al., 1996). Data are displayed on a multi-plot trellis (matrix) of graphs. A
response and an explanatory variable are on the matrix axes. Conditioning explanatory
variables are defined and each pane1 (cell) within the matrix shows the response vs. the
primary variables for a unique combination of the conditioning variables. Some features
included to enhance visual comparisons between panels are: (1) ordering g-raphs by
numerical explanatory variables; (2) standardizing axes values across panels; (3)
including grids; and (4) labeling the panels using strip labels. For this example, both
mean and individual growth curves within different sex, litter-size and treatment
categories were displayed graphically. This was first done for a11 the data and then for
individual village and herd-within-village strata to assess the consistency of the mean
growth rates and their variability within each explanatory variable category graphically
by village and herd (e.g. visually comparing the mean growth rates of lambs by sex and
litter-size class among the fïve treatment and five control villages).
R. Lancelot et ai. /Preventive Veterinary Medicine 46 (2000) 225-247
229
2.2.2. Linear-multilevel grrîwth mode1
The linear-multilevel growth mode1 is an extension of the mixed-effects models
described by Rao (1965) for growth curves and by Laird and Ware (1982) for
longitudinal-data analysis. The mode1 was specified according to the notation of Pinheiro
and Bates (2000), as follows. The highest grouping level (level 3) was lamb with
subscript k. The intermediate level (level 2) was the herd with subscript j. The lowest
grouping level (level 1) was village with subscript i.
The mode1 for the kth lamb within the jth herd within the ith village was
yijk = &jkb + Zi,&‘i + Zd,kb, + Z&‘ijk + Eijk,
i= 1,. ..,M, j - 1 >..., Mi, k=l>...> M,, biwV(O,&),
b, - N(O, &), byk N N(O, x3), eijk N N(O, dn,k),
(1)
where
l M is the number of villages: M=l0 in this study.
l Mi is the number of herds in village i: Cz,Mj= 45.
l M, is the number of lambs in the herd j in tbe village i: CzIxjTIMq = 182.
l Y$ is a vector of weights for animal k in herd j and village i, wlth the vector length
=n$ and ~E1~,~l~,,n@ = 922 live weights. This mode1 (and related fitting
algorithms) is appropriate for unbalanced data (Laird and Ware, 1982). It does not
make any restriction on the number of weights per lamb (n$), nor on the spacing
between measurements. In this analysis, standardised live-weights were used because
of data availability (not because of limitations in method or software).
l .$k is an $kxf design matrix for the fixed effects, where f is the number of fixed
effects.
o p is a vector of coefficients for the fixed effects, with length of b=$
l zijk, Zij,k, and z$ are design matrices for the village-level, herd-level and lamb-level
random effects, respectively. Their dimensions are ngkXrvillage, &jkX&,d, and
nvk x rlmb, respectively. rvillage, rherd and rli,b are the number of random effects at
the village, herd and lamb level, respectively.
l b;, b, and b+ are random vectors of coefficients for the village-level, herd-level, and
lamb-level random effects, respeCtiVely. Their 1engthS are rVillage, rherd and r]&,,
respectively.
l XijkP i- Zijkbi $ Zy,kbg i- Zijkbijk iS the linear predictor with its fixed (&j$) and
random (Zijkbi i- Zij,kbij $ Zijkbijk) parts. The fixed predictor gives the population
means, while the random predictor represents tbe deviations around these means, with
Village, herd and lamb components (Zijkbi, zQ,kb, and Zokbgk, respectively).
l cl, &, c3 and /iijk are variante-covariance matrices for the village, herd and lamb-
level random effects, and for residual error, respectively. These matrices are square and
symmetric with ranks rVillager rh&, rlamb and &$, respectively. For computational
reasons they must be positive-definite, i.e. their eigen values must be strictly positive
(Pinheiro and Bates, 2000). When needed, further constraints cari be imposed to get
more parsimonious matrices, such as diagonal or block-diagonal matrices. Each term
on the diagona1 (variante) or off-diagonal (covariance) is a random parameter,
730
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
estimated by a complex, computationally intensive and iterative process. An important
role of graphical analysis is to guide the specification of an explanatory but
parsimonious random-effects structure. If residual autocorrelation exists, it cari be
parameterized by /i,, containing a set of parameters p (Diggle et al., 1994). Again,
fitting algorithms are slow and subject to convergence problems, because a strong
numerical competition occurs during the fitting process between C- and /I-matrices
parameters (Davidian and Giltinan, 1995). Graphical analysis cari play a key role in
highlighting whether both autocorrelation and/or covariance parameters are required
and in what mode1 for-m.
l
a2 is the vector of residual variante. Additional parameters, q, to allow for non-constant
variante over the growth period cari be considered. Again, as for the estimation of
C- and /i-matrices described above, graphical analysis cari be very valuable.
Models were fitted using the restricted maximum-likelihood (REML) method
(Patterson and Thompson, 1971). Models with identical fixed effects - but different
random effects - were compared according to the Akaike’s Information’s Criterion
(Akaike, 1974), where AIC=-2 x log-likelihood+2 x npar, with npar being the number of
fixed and random parameters in the model. Mode1 likelibood (which is always higher
when adding new parameters) is penalized by the number of parameters estimated. AIC is
thus a trade-off between maximizing likelihood and minimizing the number of
parameters to be estimated. The AIC is increasingly used in regression analysis (e.g.
Venables and Ripley, 1999, pp. 185-188). Using this formulation of the AIC, the mode1
with the smallest AIC was retained.
2.2.3. Modeling strategy
Growth curves were fitted by polynomials of age, which is a common practice in
growth analysis (e.g. Rao (1965) or Lee (1991) for theory, van der Leeden et al. (1996)
for a comparison of methods and software, and Green et al. (1998), Orji and Steinbach
(1981) for applications). Fixed effects (main effects and interactions) were chosen
according to a trellis display of live-weight means against age, given sex, litter-size and
treatment. A visual inspection of the shape of the average-growth curves determined the
order for the age polynomial. Main effects and interactions were selected based on visible
between-group differences in growth-curve intercepts and slopes. Particular attention was
paid to assess differences in slopes (growth-rate differences) among the different panels,
revealing two-way (and possibly higher-order) interactions with age. These graphically
chosen fixed effects were retained in subsequent random-structure selection steps.
Trellis displays were also used to assess village-, herd- and lamb-level random effects.
As a first step, mean growth curves for each sex and litter-size category were compared
between treatment and control villages and herds within villages for consistency of effect.
Individual-live-weight trajectories were then plotted for each combination of sex, litter-
size and treatment for a11 the data and then for individual village and herd-within-village
strata. Variations in growth-curve shapes, intercepts and slopes were examined to Select
possible random effects.
A number of methods were used to assess patterns of correlation graphically over the
growing period. First, ordinary least-squares (OLS) regressions were performed on each
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
231
lamb; live-weight was the response and a polynomial of age (plus the intercept) was the
explmatory variable. Thus, q+l parameters were estimated for each lamb, where 4 was
the order of the polynomial (e.g. three parameters for a quadratic polynomial Of age: the
intercept and linear and quadratic parameters). Scatter plots of the OLS parameters were
used to assess the variability of each parameter and the between-parameter correlation.
Plots for each village and herd were also constructed and compared. These indications
were used to specify the random-effects covariance structure for the matrix Z; and
whether it might vary for different villages and herds.
Second, plots of the empirical autocorrelation function (Diggle, 1990, pp. 34-44) were
used to check the serial correlation of the residuals. Most likely correlation structures
(matrix ,4,, the residual variante) were chosen for further investigation according to these
plots.
Third, the need for an even more complex residual variante structure (allowing for
variante terms to vary over the growth period) was assessed using plots of residuals vs.
fitted values and residuals against explanatory variables.
The results of these graphical comparisons and the number of observations at each
level helped to define the strategy for comparing specific random-effects structures at
different levels. Because most observations were available at the lamb level, more-
complex random-effect structures suggested as necessary by graphical analyses could be
fitted initially and then subsequently assessed for possible simplification. For village- and
herd-level random-effect structures, fewer observations were available. Thus, parameters
that were considered to be important based on graphical analyses were tested by adding
them to nul1 models including common village and herd correlation terms. In each case,
the AIC was used to compare different random-effect structures.
Once the random-effect structure was established, the significance of the fixed effects
was assessed (with F tests) and approximate confidence intervals parameters were
obtained (using a normal approximation to the distribution of the REMI, estimators)
(Pinheiro and Bates, 2000). Plots of observed vs. fitted values were drawn to check the
goodness of fit. Normal distributional assumptions were ch=ked with quantile-quantile
plots (Cleveland, 1993).
3. Results
3.1. Mode1 structure
3.1.1. Selection offixed effects for subsequent mode@
The number of lambs in each sex, litter-size and treatment category are listed in
Table 1. Sample-mean live-weight growth curves (for each combination of the
explanatory variables: age, sex, litter-size and treatment) were plotted as a trellis display
in Fig. 1. Mean curves had simple shapes that were well fitted by linear or quadratic
(mostly for lambs bom in single litters) polynomials. A quadratic polynomial was
retained to mode1 the age-live-weight relationship in subsequent regression analysis. Sex
of lambs was not associated with growth, in any litter-size or treatment category. A litter-
size effect was observed in each treatment group, both for intercepts (main effect) and
232
R. Luncelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
Table 1
Number of kdmbs in each of the main cross-classification StI”dki (182 lambs bom in 45 herds and 10 VihgeS in
Kolda - Senegal, between October and December, 1987)
Sex
Litter-si&
Treatrnen t
Number of lambs
Female
Single
Co&ol
12
Morantel
29
Multiple
Control
18
Morantel
26
Male
Single
C o n t r o l
22
Mora&l
41
Multiple
Control
1.5
Morantel
19
a Sixty-eight twin and 10 triplet lambs were combined in the “multiple” category.
slopes (agexlitter-size interaction). For single lambs, mean growth curves were similar
for lambs in the control and morantel groups. However, for multiple lambs, the slope was
higher for the morantel than for the control groups. Because the intercepts (live-weight
at 15day of age) were close and the slopes differed, there was evidence for an agex
litter-sizextreatment interaction. We decided to retain age, litter-size and treatment
main effects and age xlitter-size, litter-size xtreatment and age xlitter-size xtreatment
Morantel
X-X-X
Corltrol
o - o - o
15 45 75
6
4
15 4s 75
&e(dws)
Fig. 1. Mean live-weight growth curves for 182 Djallonke lambs bore in 10 villages in Kolda (Senegal) between
October and December, 1987. The strip labels at the top of each pane1 indicate the sex (male of female) and the
litter-size (single or multiple) categories.
R. Lmcelot et al. /Preventive Veterinaty Medicine 46 (2000) 225-247
233
15 45 75
*
I. 10 a te: ::
15 45 75 1
Age W
Fig. 2. Individu4 live-weight gmwth curves for 182 Djallonke lambs born in 10 villages in Kolda (Senegal)
bctween October and December, 1987. The strip labels at the top of each pane1 indicate the sex (male of female),
the litter-size (single or multiple) and the treatment (contml or Morantel) categories.
interactions as the fixed effects. The agextreatment interaction was also included because
of the presence of the second-order interaction (agexlitter-sizextreatment).
3.1.2. Selection of random eflects for subsequent mode&
Individual growth curves were plotted on a trellis display (Fig. 2). Most of the
individual curves had simple pattems. A large variation was observed for the 15-day live-
weight (indicating that a random intercept should be included). Few individual growth
curves crossed each other (showing that small lambs tended to stay small and large lambs
large). This observation and changes in variante with age suggested that random slopes
should be included (as well as covariance parameters between intercept and slopes).
A trellis display of herd-averaged growth curves (not shown) revealed similar pattems
as in Fig. 2 (supporting further investigation of herd-level random intercepts, slopes and
covariance parameters). Within-village herd-averaged growth curves were plotted as a
trellis display (Fig. 3). This graph shows that while patterns of growth across villages
were roughly similar, most villages had very few herds. Thus, it was difficult to
distinguish village from herd variability and village and herd-within-village analyses
would lack power (a common difflculty in multilevel studies).
The scatter plots of OLS parameters (Fig. 4) showed large variations of the parameters
on their own scale. A strong positive correlation was found between the intercept and
linear age. Negative correlations were found between linear and quadratic age as well as
between intercept and quadratic age. For these reasons, C3 (the lamb-level random-effects
R. Lnncelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
39
10 14
”
6
bfl
-,,
:‘.
:::
:::
::
F
y 14
;z 10
6
2 15
,54575cirll;r:
I
j A:
j
,c;’
j
j LT
.:, ’
11
I
I
-1
I
,
“I
,
,
,
,
45 75
15 45 7 5
154575
15 45 75
4s 09
Fig. 3. Within-village, herd-avemged growth curves for 182 Djallonke lambs bom in 10 villages in Kolda
(Senegal) between October and December, 1987. The strip labels at the top of each pane1 indicate the village
(e.g. DEM) and the litter-size (single or multiple) categories. The Upper-case letter in the top, left corner of each
pane1 stood for the treatment category: C, contml; M, morantel.
1-Intercept
0.5 1.5
2.5
Kg lage
Fig. 4. Scattcrplot display of correlation between ordinary-least-squares
parameters. An OIS regression was ‘iïtted
for each lamb, with live-weight as the responsc and a quadrdtic polynomial for age as the explanatory variable (plus
intercept). The resulting parameters (intercepf linearage and age2) were used for the scattcrplot. The unit forage was
a 15-day period. A regression line and the corresponding p value were added to each plot to improve visual
perception of the correlation. One large outlier (pammeter for age2=-0.6) was removed before plotting.
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
235
variance-covariance matrix) was allowed to follow an unstructured pattern (with
variantes for the intercept, age and squared-age on the diagonal, and covariances
elsewhere (six parameters in total: three variantes and three covariances)).
Because of convergence problems, we started with a rather simple random structure.
An initial mode1 (fixed effects as described in Section 3.1 .l, Cs as described above and
random-effect intercept-only terms for village (Ci) and herd (X2)) had an AIC value of
1215.26. In agreement with the graphical assessment (Fig. 3), removing the village-level
random effect reduced the AIC (AIC=1213.26)
very little, indicating that a village-level
random effect explained very little of the variability in growth. A blocked-diagonal
matrix Cz was selected, with a covariance between intercept and age and no covariance
between age and squared-age, nor between intercept and squared-age. This herd-level
random structure provided a better fit than random-effect intercept-only, or unstructured
(intercept, linear and squared-age) structures (AIC=1188.16 vs. 1213.26 and 1190.40,
respectively). The lamb-level matrix C3 was left unchanged.
Residuals autocorrelation was observed at lags 1 (15 days apart), 2 (30 days apart) and
4 (60 days apart) (Fig. 5, left plot). The Perfect cor-relation at lag 0 is trivial but is
conventionally displayed to give a better visual perception of the Y-axis scale. After trying
a variety of possible serial cor-relation structures, an order-1 autoregressive [AR(l)]
cor-relation structure was retained (AIC=1141.01)
on the basis of AK improvement and
better autocorrelation pattern (in some cases, AIC was smaller than with AR(l), but plots
of the autocorrelation function revealed disastrous cor-relation patterns). The AR(l)
parameter was an estimate of a1 in the following equation:
where rijk,r and rijk, +i were the residuals at t and t-l, respectively, and E, was a white-
noise process with zero mean and fmite variante oc (Diggle, 1990). Some cor-relation
remained in the residuals at lag 1 (Fig. 5, right plot).
An inspection of the variation of residuals as a function of age (not presented) showed
that the variante of residuals was smaller at 90 days than at younger ages. Thus, a dummy
No serial correlation
AR (1) serial correlation
-
1.0
1.0
fz
‘S
0.5
0.5
ca
-G
-:l__ï-:‘-r:---
I
,
1 -0.5
01 2 3 45
F_.
5
0.0 J _______._________._._.~--~-~-~--~
_._.---.
I
0.0
------_-..._..___________
2
t-----~
2 -0.5
01 2 345
Lag (15- d periods)
Fig. 5. Plots of the autocorrclation
function for no correlation and AR(I) correlation structures for the linear-
multilevel live-weight growlh mode1 of 182 Djallonke Iambs bom in 10 villages in Kolda (Senegal) between
October and December, 1987. The dotted lines showed the 95% confidence band for the autocorrelation.
Bars
uutside the band indicated a signifïcant autocorrelation for the rcsiduals at this lag.
236
R. Lmcelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
3 1
2
3
4
5
2 3 4 5 6 7
2
4 6
8 10
2 4 6 8 10 12
4 6 8 10 12 14
4 6 8 10 12 14 16
Observed live weight (kg)
Fig. 6. Observed vs. fitted values for the linear-multilevel live-weight growth mode1 of 182 Djallonke lümbs
bom in 10 villages in Kolda (Senegal) between October and December, 1987. The scales were different for each
panel. Points on the bisecting line occur when the mode1 predicts the weight observed.
variable was created, taking the values of 1 for age=90 days, and 0 elsewhere. This
variable defined two strata for which separate residual standard errors were estimated.
The fit of the resulting mode1 was slightly better than a mode1 with a single residual error
(homoscedastic model) (AIC: 1140.28 vs. 1141.01, respectively).
3.2. Mode1 diugnostics
The goodness-of-fit of the mode1 (fitted live-weights vs. observed live-weights at each
weighing age) is presented in Fig. 6. Because live-weight showed a large variation with
age, the scales were allowed to vary from pane1 to panel, to get a better visual perception
of the correlation. As expected, the fit was better at 90 days than at any other age.
However, no systematic deviation from the bisecting line was observed and the overall
goodness-of-fit appeared good. Quantilequantile plots (not shown) demonstrated no
unusual residual pattern that might have caused us to reject the normal-distribution
assumptions.
3.3. Numericul results
3.3.1. Fixed effects
Fixed-effect parameters estimates are listed in Table 2 and the population-averaged
growth curves displayed in Fig. 7. A comparison of these results to those found in other
West Aftican studies for 90 days live-weight in Djallonke lambs is presented in Table 3.
Table 2
Fixed effects estimates and confidence intervals for the linear-multilevel growth mode1 of 182 Djallonke lambs
bom in 10 villages in Kolda (Senegal) between October and December, 1987
Fixed effect
Estimate
95% confidence inter&
p (F test)
Lower limit
Upper limit
Intercept
3.70
3.43
3.97
p<10-3
Ageb
1.64
1.45
1.83
p<10-3
Age2b
-0.09
-0.12
-0.06
p<10-3
Litter-sizec
-1.10
-1.45
-0.76
p<lO@
Treatmentd
-0.13
-0.47
0.22
0.47
Age x Iitter-size
-0.74
-0.97
-0.51
p<10-3
Age* x litter-size
0.06
0.02
0.09
p<10-3
Age x treatment
-0.02
-0.25
0.21
0.33
Age* x treatment
0.03
-0.01
0.06
0.33
Litter-size X treatment
0.46
0.03
0.90
0.04
Age x litter-size xtreatment
0.34
0.05
0.63
0.03
Age* x litter-size x treatment
-0.05
-0.09
-0.01
0.03
a Normal approximation to the distribution of the REML estimators.
b Age and age2 were the linear and quadmtic terms for age, respectively. Age was expressed in 15day unit
and the origin was set to 15 days of age.
’ Two categories: single (reference) and multiple, the latter combining twin and triple births.
d Two categories: control (reference) and morantel.
Morantel
m-m-m
Control
o-o-o
15
30 45 60
75 90
10
Age (days)
Fig. 7. Population-averagcd live-weight growth curves açcording to litter-size and treatmcnt categories,
prçdicted from the fïxed part of the lincar-multilevcl live-weight growth mode1 of 182 Djallonke lambs bom in
10 villages in Kolda (Senegal) between October and December, 1987. The strip labels ut the top of each pane1
indicated the litter-size (single or multiple) category.
23x
R. Luncelot et ul. /Preventive Vherinury Mrdicine 46 (2000) 225-247
Table 3
Compaison of 90-day live-weights observed in different studies of Djallanke larnbs growth in West A~I%X
Keference
Country
Method
Sample size Mean
Birth type
Single
Multiple
This study (control lambs)
Senegal
LMEa
182
8.2
9.6
6.3
Agyemang et al. (1991)
T h e Gambia SMb
NA’
6.8
N A
N A
Armbruster et al. (1991)
C ô t e d’ivoire L M E
359
8.4
9.1
7.7
Charray (1986)
C ô t e d’ivoire SM
1024
9.8 NA
N A
Fall et al. ( 1983)d
Senegal
LME
360
8.1
8.7
6.5
Poivey et al. (1982)
Côte d’ivoire
OLY
293 8.6
10.0
7.1
Sumberg and Mack (1985)
Nigeria
SM
173
8.8
N A
N A
Symoens and Hxdouin (1988)
Cameroon SM
N A
11.8
13.1
9.5
Tillard (1991)
Senegal
OU
209
9.5
N A
N A
a Linear mixedeffects model.
b Sample mean.
’ Not available.
d Spline interpolation from published results.
’ Ordinxy least-squares.
Al1 the fixed effects considered important in graphical analyses were signifïcant in the
final mode1 (CI=O.~~) and those considered not important graphically (sex and its
interactions with the other explanatory variables) were not significant.
A comparison of the intercepts and slopes of the growth curves of the four important
treatment and litter-size categories (single vs. multiple litterxmorantel vs. control)
(Table 2 and Fig. 7) showed that the higher the intercept, the higher the parameter for age,
and the lower the parameter for squared-age. This latter coeffkient could be considered to
be a growth deceleration parameter.
Average daily weight gains (ADWG) in grams per day were calculated from the fixed-
effect estimates (Table 4). Initia1 differences in ADWG in the [15-30 days) age class were
associated with parameters related to linear age (agexlitter-size, agextreatment and
agexlitter-sizextreatment, see Table 2). Growth rate decreased with age. The higher the
Table 4
Averdge daily weight gains according to litter-six-by-treatment
categories, predicted by the lincar-multilevel
growth mode1 of 182 Djallonke lambs bom in 10 villages in Kolda (Senegal) between October and December,
1987 (population-averdged values)
Age ChS (dayS)
Average daiiy weight gain (grams per day)
Single lambs
Multiple lambs
Control
Morantel
Co&ol
Morantel
1 U-30)
103
103
58
78
13@-45)
91
95
53
70
[45-60)
79
86
49
63
[60-75)
67
77
45
5.5
[ 75-90)
55
68
40
48
R. Lancelot et al. /Preventive kterinary Medicine 46 (2000) 225-247
239
initial ADWG, the greater tbe decrease (i.e. growth rate decreased faster for animals with
high initial weight). These different growth decelerations were associated with
parameters related to squared-age (age x litter-size, age x treatment and age x litter-
sizex treatment, see Table 2). However, the initial treatmentxlitter-size group live-
weight rankings remained the same over a11 ages: single lambs had higher growth rates
than multiple lambs and for multiple lambs, growtb rate was higher in the morantel than
in the control group.
3.3.2. Random parameters
Random parameter estimates are listed in Table 5. The correlation was r=0.93 between
intercept weights and weights at older linear ages at the lamb level. The correlation
between the intercept and square&age weights was ~-0.60 and between the age and
squared-age weights was ~-0.26. A high positive correlation also was observed
between the intercept and linear-age weights at the herd level (~0.92).
A large overall random variation was found at the lamb-level (Table 6). The
interquartile range for the random deviation from the population mean (fixed-effect
prediction) was 2.6 kg at 90 days, and the range was 9.2 kg. At this age, population
Table 5
Random pammeters estimates for the linear-multilevel growth mode1 of 182 Djallonke lambs bom in 10 villages
in Kolda (Senegal) between October and December, 1987
Random effect
Estimate
95% confidence intervdla
Lower limit
Upper limit
Herd lrvel
S.D.(intercept)b
0.30
0.17
0.50
S.D.(age)’
0.21
0.13
0.34
Correlation (intercept, age)
0.92
0.23
0.99
S.D.(age’)”
0.04
0.02
0.05
Lamh Ievel
S.D.(intercept)
0.55
0.46
0.66
S.D.(agc)
0.33
0.27
0.40
S.D.(age2)
0.03
0.02
0.04
Correldtion(intercept,
age)
0.93
0.45
0.99
Correlation(intercept,
age*)
-0.60
-0.87
-0.06
Correlation(age, age*)
-0.26
-0.62
0.19
Order-1 nutorqressive
parameter
UI
0.6U
0.47
0.79
Hrreroscedusticiry purumeter
h
0.48
0.17
1.33
Nesidtral stnndnrd error
s
0.34
0.26
0.44
a Normal approximation to the distribution of the REML estimators.
b Standard deviation of the random effect.
’ Age was expressed in 15day units and the origin was set to 15 days of age. Age and age2 were the linear
2nd quadratic random effects for age, respectively.
240
R. Lancelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
Table 6
Live-weight differences due to the rdndom effects in the linear-multilevel gmwth mode1 of 182 Djallonke lambs
born in 10 villages in Kolda (Senegal) between October and December, 1987
Age (days)
Rdndom deviation (kg)
Herd effect
Lamb effect
Overall deviation
IQRa
Rangeb
IQR
Range
IQR
Range
1 5
0.3
0.9
0.7
2.5
0.9
2.9
30
0.5
1.5
1.1
4.2
1.4
4.8
4s
0.7
2.2
1.5
5.8
1.9
6.7
60
0.9
2.9
1.8
7.0
2.3
7.8
7.5
1.2
3.7
2.1
7.8
2.8
9.7
90
1.5
4.7
2.6
9.2
3.3
11.6
a Interquartile range was the difference (due to the random effects) between the thiid and first quartiles of
live-weight devidtions from the population mean.
b Difference (due to the random effects) between the maximum and minimum values of live-weight
deviations from the population mean.
2
0
-2
2
0
-2
1.5 52 90
15 52 90
15 52 90
15 5 2 9 0
15 5290
Age (days)
Fig. 8. Deviation from population mean due to the herd-level random effects in the linear-multilevel live-weight
growth mode1 of 182 Djallonke lambs bom in 10 villages in Kokid (Senegal) between Octobcr and December,
1987. The strip label at the top of each pane1 was tbe herd identifier. Herds were ordered from the lowest (bottom
Icft panel) to the highest (top right panel) 90-day dcviation. The Upper-case letter in the top, left corner of each
pane1 stood for the treatment category: C, control; M, morantel.
R. bncelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
2 4 1
means, predicted from the fixed part of the model, ranged from 6.3 kg (multiple lambs,
control group) to 10.0 kg (single larnbs, morantel group). The herd-level variation around
the population growth curve was also large (Table 6 and Fig. 8). Treatment category was
not associated with the 90 days ranking of herd-level random deviations (Wilcoxon test:
W=648, ~28, m=17, p=O.93).
The AR(l) parameter (ai=0.66, with age expressed in 15-day periods, Table 5)
indicated that residuals’ autocorrelation persisted over at least 1 month (two 15-day
periods). The additional heteroscedasticity parameter in Table 5 (h=0.48) indicated that
the standard error of weight estimates at 90 days was about half as large as at other ages.
4. Discussion
The analysis of complex, multilevel, repeated-measures data requires that the structure
of fixed and random-effects in the mode1 be specified appropriately for correct inferences
to be made. In this example, the main study objective was to compare the growth rates of
lambs of morantel-treated vs. untreated ewes. However, a correct inference of the effect
of morantel treatment depended on: (1) identifying important fixed effects and their
interactions that could influence the morantel-treatment effect; (2) specifying the correct
random-effects structure at the multiple levels of village, herd, lamb and measurements
within lamb; and (3) specifying the correct cor-relation structure between repeated weight
measurements for individual lambs. This paper has demonstrated that graphical methods
cari play a useful role in appropriately identifying the important fixed effects and their
interactions and in specifying a logical random-effect structure to guide complex mode1
building.
In this paper, trellis graphical methods clearly and simply identified the important fixed
effects and their interactions. Al1 frxed effects highlighted by graphical analysis were
statistically significant in subsequent models and those not highlighted were insignificant.
This was the case for the lamb-sex effect and its interactions. Without graphical support
for its exclusion, we would probably have retained it as a fixed effect with interactions,
greatly increasing the rank of the fixed-effects design matrix (X) (from 12 to 24 columns).
This is important for computing time, even with simpler OLS models. With larger data
sets, calculations with 8-10 main effects and a11 their possible interactions might
overwhelm present microcomputers.
Most importantly, simple graphs highlighted that a morantel-treatment effect only
occurred in multiple-birth lambs. This lead to a mode1 formulation that identifïed a
significant morantel-treatment effect that a previous analysis, using cross-sectional
analysis-of-variante (ANOVA), had missed (Tillard, 1991). At the time of this first
analysis, no software package was widely available to fit multilevel models. Because
of the lack of power of the cross-sectional ANOVA approach, litter-size-by-treatment,
and age-by-litter-size-by-treatment interactions were not found to be significant
CJDO.05).
Graphical methods also proved very useful at identifying random-effects structures or
potential problems with their estimation at multiple levels. This is crucial when analysts
are confronted to new kinds of data without preliminary information on correlations
242
R. Lancelot ef al. /Preventive krerinary Medicine 46 (2000) 225-247
between random effects. From our experience, it is not straightforward to choose an
appropriate (explanatory and parsimonious) structure for the random-effects variance-
covariance matrix. Cubic (or higher-order) polynomials might be necessary to mode1
live-weight growth over a longer period of time (e.g. from birth to 1 year). The number
of variance-covariance parameters in an unstructured matrix, 2, is large (10 for each
hierarchical level with a cubic polynomial) and convergence problems are likely to
occut In such situations, graphical analysis might be the only reasonable (or even
possible) way to specify a plausible and numerically stable random structure. In this
analysis, initial graphical assessments were consistently confïrmed in subsequent models.
At the herd level, growth pattems varied and this Was reflected in significant herd
random intercepts and slopes. At the lamb level, trellis plots revealed more-complex
correlation patterns (leading to model-building strategies starting with a complex
correlation structure). For individual weight measurements within lamb, graphical
methods to identify important correlations at different lags and differences in residual
variation at different ages also provided useful information. In our opinion, the graphical
analyses used were a powerful tool in guiding mode1 building and in effectively
displaying important results.
Dam might have been considered as a grouping level in this analysis. As a matter of
fact, lambs bom from multiple litters might share common and unobserved genetic,
nutritional and health features, inducing cor-relation among measurements from the same
litter. This level was tested (results not shown here) but no significant random effect nor
covariance parameter could be identifïed. This is probably because most of births were
single or of small size (182 lambs from 147 ewes and 147 litters: 104 single, 68 twin and
10 triple lambs) and/or the litter-size fixed effect accounted for within-litter correlation.
With more-prolific species such as West-African-Dwarf goats (for which we have similar
data sets) over a longer observation period (with several litters per dam), a dam effect
might be more important.
The 90-day live-weight gain results in this study were within the range of observations
made in other studies in the region (Table 3). The lowest weight gains were observed in
Gambia (Agyemang et al., 1991); however, few details were given of the study methods.
The largest weight gains were observed in northwest Cameroon (Symoens and Hardouin,
1988) in an experimental herd under improved management (housing, supplementary
feeding, herd health). Thus, their results were not likely to be representative of the usual
growth performance in the local farming system. The other estimates in Table 3 were
similar to those reported in this study - particularly the magnitude of the multiple birth
(litter-size) effect.
The effect of litter-size on weight gain was large. Multiple larnbs were 1.1 kg lighter
than single lambs at 15 days of age and their growth rate was lower (regardless of age and
treatment group). Low growth rate for multiple lambs indicated that their nutritional
needs were not met by the ewes’ milk yields. At the end of the study period, single lambs
bom from morantel-treated ewes had a slightly greater growth rate than lambs bom from
control ewes (Fig. 7). However, the important effect of morantel treatment was in
multiple-birth lambs. At 90-day of age, multiple lambs bom from treated ewes were
1.3 kg heavier, on average, than those bom from control ewes (7.6 vs. 6.3 kg,
respectively, i.e. 22% more). The morantel effect was thought to be due to a better
R. Luncelot et al:/Preventive kterinary Medicine 46 (2000) 225-247
243
milk yield for the dewormed ewes - crucial for ewes with multiple lambs (Theriez,
1984). This study took place just after the rainy season, when nutritional and health
conditions were at their best with abundant forage resources and trop residues and low
parasitic pressure. Morantel treatment might have an even greater effect in lambs born
during the second lambing peak from March to May, when ewe nutritional stress would
be greatest. Because multiple births were frequent in Kolda (43% of the lambings in this
sample), the use of morantel for deworming ewes during the rainy season should have an
important effect on flock productivity. A thorough assessment of the benefit of
deworming should take into account its effect on other. production parameters (e.g.
mortality, fertility, prolificacy) as well as the cost of treatment. This was done in a
companion study (Lesnoff et al., 2000). Using a seasonal population-dynamics modeling
approach (steady-state model), a positive effect of deworming was found, with the
financial benefit-cost ratio estimated at 3.7 11.9, 5.41 (95% confidence interval in
brackets).
However, the persistency of this positive morantel-effect on growth at older ages needs
to be assessed in further fïeld studies. Published data are available with other drugs -but
to our knowledge, not with morantel for Djallonke sheep.
At the herd and lamb levels, intercept and linear-age random effects were highly
correlated. The persistence of this cor-relation at older ages should be checked to
determine whether early weights are as good predictors for mature weight as for 90-day
weights. The major consequence would be for the selection of females. In an experiment
with Scottish Blackface sheep, Gunn (1977) demonstrated that ewe mature size, sexual
precocity and lifetime productivity cari be affected by insuffrcient early growth. However,
in Nigeria, an analysis of station records of Yankasa ewes showed no association between
sexual precocity and litter-size, birth weight and weaning weight (Osuhor et al., 1997).
However, these are not strong data (observational rather than experimental and from a
specific station).
The identification of herd and lamb random effects provide important targets for
specific studies to improve lamb growth. The presence of herd random effects is
potentially related to differences in management, feeding and herd-health practices.
T~US, studies contrasting these features in flocks with low and high growth rates
should be considered. Lamb random effects might be associated with genetic
differences, ewe nutritional and health status, lamb diseases and other factors at the
lamb level.
Graphical methods are useful in highlighting the best- and worst-performing flocks or
herds. A quick examination of Fig. 8 shows that 20% of herds (nine out of 45, bottom
line) were responsible for rnost low weights and approximately 20% for most high
weights (top line). Statistical methods that ignore these herd differences are likely to
rnake incorrect inferences on fixed effects - particularly with respect to too small
variante estimates (McDermott et ai., 1994). There are particular benefrts in better
specifying the fixed and random-effects structures with graphical methods. Random-
effect deviations from population averages are intuitively appealing in such multilevel
field studies - and if correctly specified, are more powerful at detecting differences than
rnarginal (e-g. Zeger and Liang, 1986) or generalized-least-squares (e-g. Box et al., 1994)
models (van der Leeden et al., 1996).
744
R. /~~nc&t er al. /Prevertlive Veteritm-y Medicitw 46 (2000) 225.247
Another important modeling issue to which graphical tools cari be applied is in the
decision of how parsimonious or complex a mode1 should be. For example, in these data,
autocorrelation pattems were complex. Even an unstructured L$ matrix (one parameter
for each lag pair, for a total of 21 correlation parameters) could not remove residual
autocorrelation. Graphical assessment (Fig. 2) showed that a small proportion of lambs
(less than 10%) grew irregularly and poorly and that their growth was not easily modeled
by age and age-squared random effects. Deleting these poor-growing lambs to improve
mode1 fit is not an option. Heroic measures to better mode1 their growth is probably not
wise either. The use of a simple AR(l) mode1 to parsimoniously explain residual
autocorrelation for most of the lambs, the recognition that poor-growing lambs bave hard-
to-mode1 growth curves (for which specific studies could be initiated if considered
important) and the use of graphical methods to investigate when poor mode1 fit may
influence inferences about fixed and random effects seems a sensible compromise.
Structured graphical analysis with trellis graphs and other recently developed methods are
important tools in efforts to effectively summarize and test biological processes with
sufficient and not excessive complexity.
5. Conclusion
The graphical methods presented here are an important adjunct to the consideration of
linear-multilevel growth models. Users cari specify relevant starting points for the fixed
and the random effects, check mode1 goodness of fit and assess distributional assumptions
within the same computing and graphie environment. Until recently, this process was
cumbersome and required different software packages to perform different steps
(Mazumdar et al., 1999). The use of graphical methods in mode1 building and assessment
should greatly enhance the ability of veterinary epidemiologists to explore different
random structures and to effectively communicate their fïndings.
In this example, graphical analyses greatly aided in highlighting important fixed and
random effects. Of particular importance was the positive effect of ewe deworming on the
pre-weaning live-weight growth of multiple-birth lambs. In addition, the relationships
between early weight measures and subsequent growth and the large variability of growth
between herds were nicely demonstrated by graphical anal,ysis of correlation patterns and
herd random effects. In future, we foresee that graphical tools Will become more user-
friendly, allowing for interactive outputs to explore specifrc questions of both researchers
and their clients. This should allow for more their broader application and the better
integration of results from individual studies into a more general livestock systems
context.
Acknowledgements
This work was funded by a joint ISRA-CIRAD research programme f’uthofogie er
p-oductivite’ des petits ruminants dans les systèmes traditionnels (Small-ruminants
production and health in traditional farming systems). We thank Olivier and Brigitte
R. -celot et af. /PreVentiVe Vetetinary Medicine 46 (2000) X5-247
245
Faugère for providing us with data, and ISRA, Senegalese Veterinary Services and
farmers for their co-Op=tion. We also thank José Pinheiro (Lucent Technologies,
Murray HiIl, NJ) and Douglas Bates (University of Wisconsin, Madison, WI) for their
help in handling and understanding linear-multilevel models, and for providing us wit,h an
early version of their book. The two anonymous reviewers provided helpful comments on
the fiist draft of this paper.
References
Abassa, K.I?, Pessinaba, J., Adeshola-Ishola,
A., 1992. Croissance pi-sevrage des agneaux Djallonké au centre
de Kolokopé (Togo). Rev. Élev. Méd. Vét. Pays Trop. 45 (l), 49-54.
Agyemdng, K., Rawlings, P., Clifford, D., Bojang, N., Samba, A., 1991. Productivity and health parameters of
mal1 ruminants in villages of the Gambia. Bull. Anim. Hlth. Prod. Afr. 39, 129-135.
Akaike, H., 1974. A new look at statistical mode1 identification. IEEE Trans. Automatic Control AU-19, 716-
722.
-bruster, T., Peters, K.J., Metz, T., 1991. Sheep production in the humid zone of West Africa. II. Growth
performance and liveweight of sheep in improved and h-aditional production systems in Côte d’ivoire. J.
Anim. Breed. Genet. 108, 210-219.
Becker, R.A., Cleveland, W.S., Shyu, M.-J., 1996. The visual design und control of trellis display. J. Comput.
hphicti Statist. 5 (2), 123-155.
Box, G.E.P., Jenkins, G.M., Reinsel, G.C., 1994. Time Series Analysis: Forecasting and Contml, 3rd Edition.
Prentice-Hall, Englewood Cliffs, NJ, 598 pp.
Charray, J., 1986. Performances des brebis Naines de I’Afrique de I’Ouest entretenues suivant deux rythmes
différents d’accélération de h reproduction, Rev. Élev. Méd. Vét. Pays. Trop., 39 (1) 15 l-160.
Clément, V., Poivey, J.P., Faugère, O., Tillard, E., Lancelot, R., Gueye, A., Richard, D., Bibé, B., 1997. Etude de
la variabilité des caractères de reproduction chez les petits ruminants en milieu d’élevage traditionnel au
Sénégal. Rev. Élev. Méd. Vét. Pays Trop. 50 (3), 235-249.
Cleveland, W., 1993. Visualizing Data. Summit. Hobart Press, 360 pp.
Cleveland, W., 1994. The Elements of Graphing Data. Summit. Hobart Press, 297 pp.
Davidian, M.M., Giltinan, DM., 1995. Nonlinear Models for Repeated Measurement Data. Chapman & Hall,
London, 359 pp.
Diggle, P.J., 1990. Tiie Series. A Biostatistical Introduction. Clarendon Press, New York, 257 pp.
Diggle, P.J., Liang, K.-Y., Zeger, S.L., 1994. Analysis of Longitudinal Data. Clarendon Press, Oxford,
253 pp.
Fall, A., Diop, M., Sandford, J., Gueye, E., Wissocq, Y.J., Durkin, J., Trail, J.C.M., 1983. Étude sur la
productivité moutons Djallonké au Centre de Recherches Zootechniques de Kolda, au Sénégal. 2. Poids
corporel, productivité des brebis et du troupeau. Rev. Élev. Méd. Vét. Pays Trop. 36 (3), 283-289.
Paugère, O., Faugère, B., 1986. Suivi de troupeaux et contrôle des performances individuelles des petits
ruminants en milieu traditionnel africain. Aspects méthodologiques. Rev. Élev. Méd. Vét. Pays Trop. 39 (l),
29-40.
Wugère, O., Faugère, B., 1993. Panurge: suivi individuel dans les systèmes d’élevage traditionnel. Maisons-
Alfort (Frdnce), CIRAD-EMVT, 339 pp.
Faugère, O., Dockès, A.-C., Perret, C., Faugère, B., 1990. L’élevage traditionnel des petits ruminants au Sénégal.
1. Pratiques de conduite et d’exploitation des animaux chez les éleveurs de la région de Kolda. Rev. Élev.
Méd. Vét. Pays Trop. 43 (2), 249-259.
Faugère, O., Tillard, E., Faugère, B., 1990. Prophylaxie chez les petits ruminants au Sénégal: régionalisation
d’une politique nationale de protection sanitaire. In: Rey, B., Lebbie, S., Reynolds, L. (Eds.), First Biennial
Conference of the Africun SKndll Ruminant Research Network, Nairobi, INCA, pp. 307-314.
Fri&che, T., Kaufmann, J., Mister, K., 1993. Parasite spectrum and seasonal epidemiology of gastrointestinal
nematodes of small ruminants in the Gambia. Vet. Parasitol. 49, 271-283.
R. hcelrzr et al./Preventive Veterinaty Medicine 46 (2000) 225-247
Gol&t&,, H., 1986. Multihd mixed linCX mode1 andySiS USing iterdtive generdlized least squares. Biomet&+
73 (1). 43-56.
Goldst&, H., 1995. MultiIevel Statistical Models. Arnold, London, 178 pp.
cEen, L-E., ~~~~~ E., Morgan, K.L., 1998. A multi-level mode1 of data with repeated measures of the effect
of bb d&&fja on weight. Prev. Vet. Med. 36 (2), 85-94.
Gc~, Y-T., McDermott, J.J., Schukken, Y.H., Hertl, J.A., Eicker, S.W., 1999. Analysis of con-ekited continuous
reps observations: modelling the effect of ketosis on milk yield in dairy cows. Prev. Vet. Med. 39, I37-
153.
Gunn, R.C., 1977. The effects of two nutritional envimnments from 6 weeks pre partum to 12 months of age on
lifetime performance and reproductive potential of Scottish Bhu&face ewes in two adult environments.
Anim. Prod. 25, 1.55-164.
Kraft, I., de L.eeuw, J., van der Leeden, R., 1994. Review of fïve rnultilevel analysis pmgrams : BMDP-SV,
GENMOD, HLM, h4L3, VARCL. The Am. Statist. 48 (4), 324-335.
Laird, N.M., Ware, J.H., 1982. Random-effectb: models for longitudinal data. Biometrics 38, 963-974.
Lancelot, R., Faye, B., Ju&s, X., Ndiaye, M., Pérochon, L., Tillard, E., 1998. La base de données BAOBAB: un
outil pour modéliser la production et la santé des petits ruminants dans les systèmes d’élevage traditionnels
au Sénégal. Rev. Élev. Méd. Vét. Pays Trop. 51 (2), 135-146.
Lee, JC., 1991. Tests and mode1 selection for the general growth curve model. Biometrics 47, 147-159.
Lesnoff, M., Lancelot, R., Tillard, E., Dohoo, LR., 2000. A steady-state appmach of benefit-cost analysis with a
periodic Leslie-m&% model: presentation and application to the evduation of a sheepdiseases preventive
scheme in Kolda Senegal. Prev. Vet. Med. 46, 113-128.
Littel, R.C., MilIiken, C.A., Stroup. W.W., Wolfimger, R.D., 1996. SAS System for Mixed Models. Cary, NC,
663 PP.
Mazumdar, S., Begley, A.E., Houck, P.R., Yang III, Y., Kupfer, D.J., 1999. Residual analysis in random
regrcssion using SAS and S-PLUS. Comput. Meth. Prograns Biomed. 58, 281-282.
McDermott, JJ., Schukken, Y.-H., Shoukri, M.M., 1994. Study design and analytic methods for data collected
from clusters of animals. Prev. Vet. Med. 18, 17.5-191.
McDermott, J.J., O’Callaghan, C., Kadohira, M., 1997. Experiences in multi-level models for assessing the
distribution of animal diseases. Épidémiol. Santé Anim. 31/32, 13.23.1-13.23.3.
Orji, B.I., Steinbach, J., 1981. Post-weaning growth and development of Nigerian dwarf sheep. Trop. Anim.
Hlth. Pmd. 13, 101-106.
Osuhor, C.U., Osinowo, O.A., Nwagu, B.I., Dennar, F.O., Abdullahi-Adee,
A., 1997. Factors affecting age at
fmt lambing in Yankasa ewes. Anim. Reprod. Sci. 47, 205-209.
Patterson, H.D., Thompson, R., 1971. Recovery of inter-block information when block sizes are unequaI.
Biometrika 58, 545-554.
Pinheiro, J.C., Bates, D.M., 1996. Unrestricted parameterizations
for variance~ovariance
matrices. Statist.
Comput. 6, 289-296.
Pinheiro, J.C., Bates, D.M., 2000. Mixed-effects Models in S and S-Plus. Springer, New York, 598 pp.
Poivey, J.P., Landais, E., Berger, Y., 1982. Etude et amélioration génétique de la croissance des agneaux
Djallonké. Résultats obtenus au Centre de Recherches Zootechniques de Bouaké (Côte d’ivoire). Rev. Élev.
Méd. Vét. Pays Trop. 35 (4), 421-433.
Puyalto, C., Sanaa, M., N’Djoya, A., Planchenault, D., 1997. Factors influencing daily weight gain on Young
zebu in Garoua district, Camemon. Épidémiol. Santé Anim. 31/32, 02.03.1-02.03.3.
Rao, C.R., 1965. The theory of least squares when the parameters are stochastic and its application to t.he
analysis of growth curve. Biometrika 52 (3/4), 447-458.
Sumberg, J.E., Mack, S.D., 1985. Village production of West African Dwxf go& ad sheep in Nigeria
Trop. Anim. Hlth. Prod. 17, 135-140.
Symoens, C., Hardouin, J., 1988. Le mouton Djallonké en élevage extensif dans le Nord-Ouest Cameroun. Rev.
Élev. Méd. Vét. Pays Trop. 41 (4), 449-458.
van der Leeden, R., Vrijburg, K., de Leeuw, J., 1996. A review of two different approaches for tbe andysis of
growth data using longitudinal mixed linear models: comparing hierarchical linear regression (ML3, HLM)
and repeated maures designs with stmctured covtiance matrices (BMDPSV), The Statist. Software
Newslett. 20, 583-605.
R. Lmcelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
247
Venables, W.N., Ripley, B.D., 1999. Modem Applied Statistics with S-Plus. Springer, New York, SO1 pp.
Theriez, M., 1984. Influence de l’alimentation sur les performances de reproduction des ovins. In: Neuvièmes
Journées de la Recherche Ovine et Caprine, Paris, INRA-ITOVIC, pp. 294-326.
Tillard, E., 1991. Évaluation technico-économique des effets des prophylaxies en milieu villageois chez les
petits ruminants de la région de Kolda (Sénégal). Thèse Méd. Vét. ENVA. Faculté de Médecine de Créteil,
172 pp.
Tuah, A.K., Baah, J., 1985. Repmductive performance, pre-weaning growth rate and pre-weaning lamb
mortality of Djallonke sheep in Ghana. Trop. Anim. Hlth. Prod. 17, 107-113.
Zeger, S.L., Liang, K.Y., 1986. Longitudinal data analysis for discrete and continuous outcomes. Biometrics 42,
121-130.
/
i/
PREVENTIVE
VETERINARY
MFD1cINE
ET ‘?F;VTI;R
Pwventive Velerinarv
I_IU”U
. .&.a*.
_
.-.-
_.__.-
Medicine 46
__......-,
(2000) 225-247www.ekevier.nl/locate/prevetmed
Graphical approaches to support the analysis of
linear-multilevel models of latib pre-weaning
growth in Kolda (Senegal)
Renaud Lancelota7*, Matthieu LesnoffQ, Emmanuel Tillardb,
John J. McDerrnottc
%Zentre de Coopération Internationale en Recherche Agmnomique pour le Développement, Département
Élevage et Médecine Vétérinaire (CIRAD-EMVT),
TA 3O/A, Campus International de Baillarguet, 34398
Montpellier Cedex 5, France
bCIRAD-Eh4Vr Station de Ligne-Paradis, 7 chemin de I’IRAT, Ligne-Paradis, 97410
St Pierre, Ile de la Réunion, France
cIntemationai Livestock Research Institute (ILRI), PO Box 30709, Nairobi, Kenya
Received 28 December 1999; accepted 30 May 2000
Abstract
Linear-multilevel models (LMM) are mixed-effects models in which several levels of grouping
may be specified (village, herd, animal, . . .). This study highlighted the usefulness of graphical
methocis in their analysis through: (1) the choice of the fïxed and random effects and their structure,
(2) the assessment of goodness-of-fit and (3) distributional assumptions for random effects and
residuals.
An LMM was developed to study the effect of ewe deworming with morantel on lamb pre-
weaning growth in a field experiment involving 182 lambs in 45 herds and 10 villages in Kolda,
Senegal. Growth was described as a quadratic polynomial of age. Other covariates were sex, litter-
size and treatment. The choice of fïxed and random effects relied on three graphs: (1) a trellis
display of mean live-weight vs. age, to Select main effects and interactions (fixed effects); (2) a
trellis display of individual growth curves, to decide which growth-curve terms should be included
as random effects and (3) a scatter plot of parameters of lamb-specifk regressions (live-weight vs.
quadratic polynomial of age) to choose the random-effects covariance structure.
Age, litter-size, agexlitter-size,
litter-size x treatment and age xlitter-size xtreatment were
selected graphically as fïxed effects and were signifïcant (~~0.05) in subsequent statistical
models. The selection of random-effect structures was guided by graphkal assessment and
* Corresponding author. Present address: Institut Sénégalais de Recherches Agricoles, Laboratoire National
d’Élevage et de Recherche Vétérinaire (ISRA-LNERV), BP 2057, Dakar-Hann, Sénégal. Tel.: +221-832-49-02;
fax: +221-821-18-79.
E-muil uddress: renaud.lancelot@cirad.fr (Ii. Lancelot).
0167-5877/00/$ - see front matter 0 2000 Elsevier Science B.V. AL1 rights reserved.
PII: SO167-5877(00)00155-0
226
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
comparison of the Akaike’s information criterion for different models. The final random-effects
selected included no random effect at tbe village level but intercept, age and squared-age at the
herd and lamb levels. Tbe structure of tbe random-effects variance-covariance matrices were
blocked-diagonal at the herd level and unstructured at the lamb level. An order-1 autoregressive
structure was retained to account for serial correlations of residuals. Smaller residual variante at
90 days than at younger ages was modeled witb a dummy variable taking a value of 1 at 90 days
and 0 eIsewhere.
Ewe-deworming with morantel during tbe rainy season lead to higber latnb live-weights
(probably related to a better ewe-nutrition and -health status). A positive correlation was
demonstrated between early weight and growth rate at the population level (with important lamb
and herd-level random deviations). Tbe persistence of .this correlation at older ages should be
checked to determine whether early weigbts are good predictors of mature weights and ewe-
reproductive lifetime performance. 0 2000 Elsevier Science B.V. Al1 rights reserved.
Keyworuk: Lmear-multievel growth model; Random effect; Grdphical methods; hz-weaning growth, Castro-
intestins parasitism; Mordntel; Sheep; Senegd
1. Introduction
Studies in veterinary medicine frequently measure phenomena repeated over time
(such as live-weights of the same animais). The analysis of these measures as multivariate
outcomes over time is usually more informative than repeated cross-sectjonal studies, e.g.
live-weights at given ages or cumulated production data such as 305day milk yield for
dairy cows (Grohn et al., 1999). In addition, veterinary studies ofien rely on multilevel
designs, in which a number of units at different levels (for example: villages, herds within
villages, and animals within herds) are sampled. Data arising from these studies often
have both longitudinal and multilevel features. Mixed-effect models have been used to
analyze such longitudinal and multilevel data appropriately (Diggle et al., 1994). They
provide a choice of effect type (either fixed or random) depending on study design and
objectives at different grouping levels. The statistical properties of linear-multilevel
models have been well described in recent years (e.g. Goldstein, 1986, 1995’; Pinheiro and
Bates, 1996) and different software developed to fit them (Kreft et al., 1994; Littel et al.,
1996).
Linear-multilevel growth models have been applied recently in veterinary studies (e.g.
Puyalto et al., 1997; Green et al., 1998). However, because of their complexity and the
different inferences that cari result from different constructs of fixed and random effects
(McDermott et al., 1997), tare must be taken in mode1 specification. Model-specification
issues include the choice of fixed and random effects at different levels and the structure
of random-effect variance-covariance matrices. The recent availabjlity of software
combining computing tools and graphjcal methods, provides new and intuitive tools for
better exploring the most-appropriate mode1 structures prior to fitting and in assessing
and comparing fitted models (Pinheiro and Bates, 2000).
In thjs paper, we emphasize the use of graphical methods as a tool in mode1 building
and assessment for repeated and multilevel outcomes. As an example, a linear-multilevel
R. Lancelot et al. /Preventive Veterinmy Medicine 46 (2000) 225-247
227
growth mode1 was developed to study the effect of ewe deworming with morantel on
lamb pre-weaning growth in a field experiment in several herds and villages in Kolda,
Senegal. Graphical and statistical analyses focused on the following objectives: (1) to
estimate the effect of morantel treatment on lamb growth correctly, (2) to derive
population-averaged estimates for lamb growth curves, (3) to evaluate the importance of
lamb, herd and village-level deviations around the populations values and (4) to provide
information on factors influencing live-weight growth in individual lambs under field
conditions. The advantages of including graphical methods in mode1 building and
assessment are highlighted.
2. Materials and methods
2.1. Data
Kolda is located in the sub-humid, southem zone of Senegal. From 1983 to 1995, mean
annual rainfall was 963 mm (mainly occurring between June and November). Means of
maximum- and minimum-daily temperatures were 35.5 and 20.6”C, respectively. The
landscape was a mosaic of forested plateaus and lowlands cultivated with rice, cotton,
maize, and groundnuts. A detailed description of the farming systems in this area was
given by Faugère et al. (1990). Briefly, the main ethnie group was the Fulani. Crop
farming was the dominant activity but livestock production was widely practiced (cattle,
sheep and goats). Sheep were mainly Djallonke (a trypanotolerant West African dwarf
breed reared extensively for meat). Normally, they were left to graze freely in the bush
except during the rainy season when they were tied to avoid trop damage. A high
incidence of gastrointestinal helminthoses occurred in sheep during the rainy season
(Tillard, 1991). An abattoir survey revealed that Trichostrongylus colubrifomis,
Oesophagostomum columbianum, Haemonchus contortus and Strongyloides papillosus
were the most-prevalent parasites in this region (Fritsche et al., 1993).
A field experiment was designed to assess the profitability of deworming sheep flocks
in this zone and was conducted from July 1987 to June 1988. Partial results were
published previously (Faugère et al., 1990; Tillard, 1991). Ten villages were selected
according to agro-ecosystem representativeness and accessibility criteria. Herds
belonging to volunteer farmers were enrolled and all their animals ear-tagged (Faugère
and Faugère, 1986). Trained surveyors recorded demographic events and growth data
during fortnightly visits. Al1 10 villages enrolled in the study were randomly allocated to
one of two treatment groups (morantel or control) using a random-number table. Five
villages were allocated to the morantel group, and five villages to the control group. Forty-
five herds were included in the experiment: 17 in the control villages and 28 in the
morantel villages. Treatments were applied at the village level. In the morantel-treated
villages, ail sheep older than 3 months were drenched three times during the rainy season
with an oral administration of morantel (Exhelm@, Pfizer) at 7.5 mg kg-’ live-weight.
No sheep in the control villages and no lamb younger than 3 months in the morantel
villages were given any dewormer. Thus, the treatment-effect on lambs was mediated
through its effect on their ewe.
228
R. Lancelot et a[. /Preventive Veterinaty Medicine 46 (2000) 225-247
In Kolda, lambing distribution is bimodal, with two peaks from October to December
(end of the rainy season) and March to May (end of the hot, dry season) (Clément et al.,
1997). Tbese two lambing seasons define two cohorts of lambs subject to different
nutritional and environmental conditions under different husbandry practices (Faugère
et al., 1990), all with possible important effects on live-weight growth. Because the
purpose of this paper was mainly methodological, we decided to limit the study to the
cohort of lambs born from 1 October to 31 December 1987. A more-comprehensive study
of the whole data set Will be presented in subsequent papers.
There were 182 lambs with 922 weight measurements in this data subset. At the time of
data collection (1987-1988), follow-up information was checked and stored in temporary
files (Faugère and Faugère, 1993). Weights stored in these files had been adjusted for the
standard ages of 15, 30,45, 60,715 and 90 days (raw weights were thrown away). Before
analysis, these files were gathered in a relational database (Lancelot et al., 1998). Given
that each lamb was weighed every 2 weeks, little information was lost on the curvature of
live-weight-growth trajectories. However, because of this standardization, a few weights
(6%) were discarded. We believe this loss of information had no important effect on our
results, beyond a small decrease in statistical power.
Previous studies suggested that sex and litter-size (number of lambs bom) were
important factors to consider for lamb pre-weaning growth (Poivey et al., 1982; Fall et al.,
1983; Tuah and Baah, 1985; Armbruster et al., 1991; Abassa et al., 1992). T~US, sex
(male or female), litter-size (single or multiple: twins and triplets) and treatment
(morantel or control) were the flxed-effect classification variables used.
2.2. Data analysis
2.2-l. Craphical analyses
Graphical analyses followed the three-step exploratory data-analysis (EDA) framework
proposed by Cleveland (1993): (1) display the data with respect to its inherent structure
and the mode1 to be fitted; (2) display the fitted mode1 and check its assumptions and (3)
examine mode1 consistency with the original data. TO enhance visual comparisons of
different variable combinations within hierarchically structured data, Cleveland (1994)
proposed the method of trellis graphies for which software was subsequently developed
(Becker et al., 1996). Data are displayed on a multi-plot trellis (matrix) of graphs. A
response and an explanatory variable are on the matrix axes. Conditioning explanatory
variables are defined and each pane1 (cell) within the matrix shows the response vs. the
primary variables for a unique combination of the conditioning variables. Some features
included to enhance visual comparisons between panels are: (1) ordering g-raphs by
numerical explanatory variables; (2) standardizing axes values across panels; (3)
including grids; and (4) labeling the panels using strip labels. For this example, both
mean and individual growth curves within different sex, litter-size and treatment
categories were displayed graphically. This was first done for a11 the data and then for
individual village and herd-within-village strata to assess the consistency of the mean
growth rates and their variability within each explanatory variable category graphically
by village and herd (e.g. visually comparing the mean growth rates of lambs by sex and
litter-size class among the fïve treatment and five control villages).
R. Lancelot et ai. /Preventive Veterinary Medicine 46 (2000) 225-247
229
2.2.2. Linear-multilevel grrîwth mode1
The linear-multilevel growth mode1 is an extension of the mixed-effects models
described by Rao (1965) for growth curves and by Laird and Ware (1982) for
longitudinal-data analysis. The mode1 was specified according to the notation of Pinheiro
and Bates (2000), as follows. The highest grouping level (level 3) was lamb with
subscript k. The intermediate level (level 2) was the herd with subscript j. The lowest
grouping level (level 1) was village with subscript i.
The mode1 for the kth lamb within the jth herd within the ith village was
yijk = &jkb + Zi,&‘i + Zd,kb, + Z&‘ijk + Eijk,
i= 1,. ..,M, j - 1 >..., Mi, k=l>...> M,, biwV(O,&),
b, - N(O, &), byk N N(O, x3), eijk N N(O, dn,k),
(1)
where
l M is the number of villages: M=l0 in this study.
l Mi is the number of herds in village i: Cz,Mj= 45.
l M, is the number of lambs in the herd j in tbe village i: CzIxjTIMq = 182.
l Y$ is a vector of weights for animal k in herd j and village i, wlth the vector length
=n$ and ~E1~,~l~,,n@ = 922 live weights. This mode1 (and related fitting
algorithms) is appropriate for unbalanced data (Laird and Ware, 1982). It does not
make any restriction on the number of weights per lamb (n$), nor on the spacing
between measurements. In this analysis, standardised live-weights were used because
of data availability (not because of limitations in method or software).
l .$k is an $kxf design matrix for the fixed effects, where f is the number of fixed
effects.
o p is a vector of coefficients for the fixed effects, with length of b=$
l zijk, Zij,k, and z$ are design matrices for the village-level, herd-level and lamb-level
random effects, respectively. Their dimensions are ngkXrvillage, &jkX&,d, and
nvk x rlmb, respectively. rvillage, rherd and rli,b are the number of random effects at
the village, herd and lamb level, respectively.
l b;, b, and b+ are random vectors of coefficients for the village-level, herd-level, and
lamb-level random effects, respeCtiVely. Their 1engthS are rVillage, rherd and r]&,,
respectively.
l XijkP i- Zijkbi $ Zy,kbg i- Zijkbijk iS the linear predictor with its fixed (&j$) and
random (Zijkbi i- Zij,kbij $ Zijkbijk) parts. The fixed predictor gives the population
means, while the random predictor represents tbe deviations around these means, with
Village, herd and lamb components (Zijkbi, zQ,kb, and Zokbgk, respectively).
l cl, &, c3 and /iijk are variante-covariance matrices for the village, herd and lamb-
level random effects, and for residual error, respectively. These matrices are square and
symmetric with ranks rVillager rh&, rlamb and &$, respectively. For computational
reasons they must be positive-definite, i.e. their eigen values must be strictly positive
(Pinheiro and Bates, 2000). When needed, further constraints cari be imposed to get
more parsimonious matrices, such as diagonal or block-diagonal matrices. Each term
on the diagona1 (variante) or off-diagonal (covariance) is a random parameter,
730
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
estimated by a complex, computationally intensive and iterative process. An important
role of graphical analysis is to guide the specification of an explanatory but
parsimonious random-effects structure. If residual autocorrelation exists, it cari be
parameterized by /i,, containing a set of parameters p (Diggle et al., 1994). Again,
fitting algorithms are slow and subject to convergence problems, because a strong
numerical competition occurs during the fitting process between C- and /I-matrices
parameters (Davidian and Giltinan, 1995). Graphical analysis cari play a key role in
highlighting whether both autocorrelation and/or covariance parameters are required
and in what mode1 for-m.
l
a2 is the vector of residual variante. Additional parameters, q, to allow for non-constant
variante over the growth period cari be considered. Again, as for the estimation of
C- and /i-matrices described above, graphical analysis cari be very valuable.
Models were fitted using the restricted maximum-likelihood (REML) method
(Patterson and Thompson, 1971). Models with identical fixed effects - but different
random effects - were compared according to the Akaike’s Information’s Criterion
(Akaike, 1974), where AIC=-2 x log-likelihood+2 x npar, with npar being the number of
fixed and random parameters in the model. Mode1 likelibood (which is always higher
when adding new parameters) is penalized by the number of parameters estimated. AIC is
thus a trade-off between maximizing likelihood and minimizing the number of
parameters to be estimated. The AIC is increasingly used in regression analysis (e.g.
Venables and Ripley, 1999, pp. 185-188). Using this formulation of the AIC, the mode1
with the smallest AIC was retained.
2.2.3. Modeling strategy
Growth curves were fitted by polynomials of age, which is a common practice in
growth analysis (e.g. Rao (1965) or Lee (1991) for theory, van der Leeden et al. (1996)
for a comparison of methods and software, and Green et al. (1998), Orji and Steinbach
(1981) for applications). Fixed effects (main effects and interactions) were chosen
according to a trellis display of live-weight means against age, given sex, litter-size and
treatment. A visual inspection of the shape of the average-growth curves determined the
order for the age polynomial. Main effects and interactions were selected based on visible
between-group differences in growth-curve intercepts and slopes. Particular attention was
paid to assess differences in slopes (growth-rate differences) among the different panels,
revealing two-way (and possibly higher-order) interactions with age. These graphically
chosen fixed effects were retained in subsequent random-structure selection steps.
Trellis displays were also used to assess village-, herd- and lamb-level random effects.
As a first step, mean growth curves for each sex and litter-size category were compared
between treatment and control villages and herds within villages for consistency of effect.
Individual-live-weight trajectories were then plotted for each combination of sex, litter-
size and treatment for a11 the data and then for individual village and herd-within-village
strata. Variations in growth-curve shapes, intercepts and slopes were examined to Select
possible random effects.
A number of methods were used to assess patterns of correlation graphically over the
growing period. First, ordinary least-squares (OLS) regressions were performed on each
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
231
lamb; live-weight was the response and a polynomial of age (plus the intercept) was the
explmatory variable. Thus, q+l parameters were estimated for each lamb, where 4 was
the order of the polynomial (e.g. three parameters for a quadratic polynomial Of age: the
intercept and linear and quadratic parameters). Scatter plots of the OLS parameters were
used to assess the variability of each parameter and the between-parameter correlation.
Plots for each village and herd were also constructed and compared. These indications
were used to specify the random-effects covariance structure for the matrix Z; and
whether it might vary for different villages and herds.
Second, plots of the empirical autocorrelation function (Diggle, 1990, pp. 34-44) were
used to check the serial correlation of the residuals. Most likely correlation structures
(matrix ,4,, the residual variante) were chosen for further investigation according to these
plots.
Third, the need for an even more complex residual variante structure (allowing for
variante terms to vary over the growth period) was assessed using plots of residuals vs.
fitted values and residuals against explanatory variables.
The results of these graphical comparisons and the number of observations at each
level helped to define the strategy for comparing specific random-effects structures at
different levels. Because most observations were available at the lamb level, more-
complex random-effect structures suggested as necessary by graphical analyses could be
fitted initially and then subsequently assessed for possible simplification. For village- and
herd-level random-effect structures, fewer observations were available. Thus, parameters
that were considered to be important based on graphical analyses were tested by adding
them to nul1 models including common village and herd correlation terms. In each case,
the AIC was used to compare different random-effect structures.
Once the random-effect structure was established, the significance of the fixed effects
was assessed (with F tests) and approximate confidence intervals parameters were
obtained (using a normal approximation to the distribution of the REMI, estimators)
(Pinheiro and Bates, 2000). Plots of observed vs. fitted values were drawn to check the
goodness of fit. Normal distributional assumptions were ch=ked with quantile-quantile
plots (Cleveland, 1993).
3. Results
3.1. Mode1 structure
3.1.1. Selection offixed effects for subsequent mode@
The number of lambs in each sex, litter-size and treatment category are listed in
Table 1. Sample-mean live-weight growth curves (for each combination of the
explanatory variables: age, sex, litter-size and treatment) were plotted as a trellis display
in Fig. 1. Mean curves had simple shapes that were well fitted by linear or quadratic
(mostly for lambs bom in single litters) polynomials. A quadratic polynomial was
retained to mode1 the age-live-weight relationship in subsequent regression analysis. Sex
of lambs was not associated with growth, in any litter-size or treatment category. A litter-
size effect was observed in each treatment group, both for intercepts (main effect) and
232
R. Luncelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
Table 1
Number of kdmbs in each of the main cross-classification StI”dki (182 lambs bom in 45 herds and 10 VihgeS in
Kolda - Senegal, between October and December, 1987)
Sex
Litter-si&
Treatrnen t
Number of lambs
Female
Single
Co&ol
12
Morantel
29
Multiple
Control
18
Morantel
26
Male
Single
C o n t r o l
22
Mora&l
41
Multiple
Control
1.5
Morantel
19
a Sixty-eight twin and 10 triplet lambs were combined in the “multiple” category.
slopes (agexlitter-size interaction). For single lambs, mean growth curves were similar
for lambs in the control and morantel groups. However, for multiple lambs, the slope was
higher for the morantel than for the control groups. Because the intercepts (live-weight
at 15day of age) were close and the slopes differed, there was evidence for an agex
litter-sizextreatment interaction. We decided to retain age, litter-size and treatment
main effects and age xlitter-size, litter-size xtreatment and age xlitter-size xtreatment
Morantel
X-X-X
Corltrol
o - o - o
15 45 75
6
4
15 4s 75
&e(dws)
Fig. 1. Mean live-weight growth curves for 182 Djallonke lambs bore in 10 villages in Kolda (Senegal) between
October and December, 1987. The strip labels at the top of each pane1 indicate the sex (male of female) and the
litter-size (single or multiple) categories.
R. Lmcelot et al. /Preventive Veterinaty Medicine 46 (2000) 225-247
233
15 45 75
*
I. 10 a te: ::
15 45 75 1
Age W
Fig. 2. Individu4 live-weight gmwth curves for 182 Djallonke lambs born in 10 villages in Kolda (Senegal)
bctween October and December, 1987. The strip labels at the top of each pane1 indicate the sex (male of female),
the litter-size (single or multiple) and the treatment (contml or Morantel) categories.
interactions as the fixed effects. The agextreatment interaction was also included because
of the presence of the second-order interaction (agexlitter-sizextreatment).
3.1.2. Selection of random eflects for subsequent mode&
Individual growth curves were plotted on a trellis display (Fig. 2). Most of the
individual curves had simple pattems. A large variation was observed for the 15-day live-
weight (indicating that a random intercept should be included). Few individual growth
curves crossed each other (showing that small lambs tended to stay small and large lambs
large). This observation and changes in variante with age suggested that random slopes
should be included (as well as covariance parameters between intercept and slopes).
A trellis display of herd-averaged growth curves (not shown) revealed similar pattems
as in Fig. 2 (supporting further investigation of herd-level random intercepts, slopes and
covariance parameters). Within-village herd-averaged growth curves were plotted as a
trellis display (Fig. 3). This graph shows that while patterns of growth across villages
were roughly similar, most villages had very few herds. Thus, it was difficult to
distinguish village from herd variability and village and herd-within-village analyses
would lack power (a common difflculty in multilevel studies).
The scatter plots of OLS parameters (Fig. 4) showed large variations of the parameters
on their own scale. A strong positive correlation was found between the intercept and
linear age. Negative correlations were found between linear and quadratic age as well as
between intercept and quadratic age. For these reasons, C3 (the lamb-level random-effects
R. Lnncelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
39
10 14
”
6
bfl
-,,
:‘.
:::
:::
::
F
y 14
;z 10
6
2 15
,54575cirll;r:
I
j A:
j
,c;’
j
j LT
.:, ’
11
I
I
-1
I
,
“I
,
,
,
,
45 75
15 45 7 5
154575
15 45 75
4s 09
Fig. 3. Within-village, herd-avemged growth curves for 182 Djallonke lambs bom in 10 villages in Kolda
(Senegal) between October and December, 1987. The strip labels at the top of each pane1 indicate the village
(e.g. DEM) and the litter-size (single or multiple) categories. The Upper-case letter in the top, left corner of each
pane1 stood for the treatment category: C, contml; M, morantel.
1-Intercept
0.5 1.5
2.5
Kg lage
Fig. 4. Scattcrplot display of correlation between ordinary-least-squares
parameters. An OIS regression was ‘iïtted
for each lamb, with live-weight as the responsc and a quadrdtic polynomial for age as the explanatory variable (plus
intercept). The resulting parameters (intercepf linearage and age2) were used for the scattcrplot. The unit forage was
a 15-day period. A regression line and the corresponding p value were added to each plot to improve visual
perception of the correlation. One large outlier (pammeter for age2=-0.6) was removed before plotting.
R. Lancelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
235
variance-covariance matrix) was allowed to follow an unstructured pattern (with
variantes for the intercept, age and squared-age on the diagonal, and covariances
elsewhere (six parameters in total: three variantes and three covariances)).
Because of convergence problems, we started with a rather simple random structure.
An initial mode1 (fixed effects as described in Section 3.1 .l, Cs as described above and
random-effect intercept-only terms for village (Ci) and herd (X2)) had an AIC value of
1215.26. In agreement with the graphical assessment (Fig. 3), removing the village-level
random effect reduced the AIC (AIC=1213.26)
very little, indicating that a village-level
random effect explained very little of the variability in growth. A blocked-diagonal
matrix Cz was selected, with a covariance between intercept and age and no covariance
between age and squared-age, nor between intercept and squared-age. This herd-level
random structure provided a better fit than random-effect intercept-only, or unstructured
(intercept, linear and squared-age) structures (AIC=1188.16 vs. 1213.26 and 1190.40,
respectively). The lamb-level matrix C3 was left unchanged.
Residuals autocorrelation was observed at lags 1 (15 days apart), 2 (30 days apart) and
4 (60 days apart) (Fig. 5, left plot). The Perfect cor-relation at lag 0 is trivial but is
conventionally displayed to give a better visual perception of the Y-axis scale. After trying
a variety of possible serial cor-relation structures, an order-1 autoregressive [AR(l)]
cor-relation structure was retained (AIC=1141.01)
on the basis of AK improvement and
better autocorrelation pattern (in some cases, AIC was smaller than with AR(l), but plots
of the autocorrelation function revealed disastrous cor-relation patterns). The AR(l)
parameter was an estimate of a1 in the following equation:
where rijk,r and rijk, +i were the residuals at t and t-l, respectively, and E, was a white-
noise process with zero mean and fmite variante oc (Diggle, 1990). Some cor-relation
remained in the residuals at lag 1 (Fig. 5, right plot).
An inspection of the variation of residuals as a function of age (not presented) showed
that the variante of residuals was smaller at 90 days than at younger ages. Thus, a dummy
No serial correlation
AR (1) serial correlation
-
1.0
1.0
fz
‘S
0.5
0.5
ca
-G
-:l__ï-:‘-r:---
I
,
1 -0.5
01 2 3 45
F_.
5
0.0 J _______._________._._.~--~-~-~--~
_._.---.
I
0.0
------_-..._..___________
2
t-----~
2 -0.5
01 2 345
Lag (15- d periods)
Fig. 5. Plots of the autocorrclation
function for no correlation and AR(I) correlation structures for the linear-
multilevel live-weight growlh mode1 of 182 Djallonke Iambs bom in 10 villages in Kolda (Senegal) between
October and December, 1987. The dotted lines showed the 95% confidence band for the autocorrelation.
Bars
uutside the band indicated a signifïcant autocorrelation for the rcsiduals at this lag.
236
R. Lmcelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
3 1
2
3
4
5
2 3 4 5 6 7
2
4 6
8 10
2 4 6 8 10 12
4 6 8 10 12 14
4 6 8 10 12 14 16
Observed live weight (kg)
Fig. 6. Observed vs. fitted values for the linear-multilevel live-weight growth mode1 of 182 Djallonke lümbs
bom in 10 villages in Kolda (Senegal) between October and December, 1987. The scales were different for each
panel. Points on the bisecting line occur when the mode1 predicts the weight observed.
variable was created, taking the values of 1 for age=90 days, and 0 elsewhere. This
variable defined two strata for which separate residual standard errors were estimated.
The fit of the resulting mode1 was slightly better than a mode1 with a single residual error
(homoscedastic model) (AIC: 1140.28 vs. 1141.01, respectively).
3.2. Mode1 diugnostics
The goodness-of-fit of the mode1 (fitted live-weights vs. observed live-weights at each
weighing age) is presented in Fig. 6. Because live-weight showed a large variation with
age, the scales were allowed to vary from pane1 to panel, to get a better visual perception
of the correlation. As expected, the fit was better at 90 days than at any other age.
However, no systematic deviation from the bisecting line was observed and the overall
goodness-of-fit appeared good. Quantilequantile plots (not shown) demonstrated no
unusual residual pattern that might have caused us to reject the normal-distribution
assumptions.
3.3. Numericul results
3.3.1. Fixed effects
Fixed-effect parameters estimates are listed in Table 2 and the population-averaged
growth curves displayed in Fig. 7. A comparison of these results to those found in other
West Aftican studies for 90 days live-weight in Djallonke lambs is presented in Table 3.
Table 2
Fixed effects estimates and confidence intervals for the linear-multilevel growth mode1 of 182 Djallonke lambs
bom in 10 villages in Kolda (Senegal) between October and December, 1987
Fixed effect
Estimate
95% confidence inter&
p (F test)
Lower limit
Upper limit
Intercept
3.70
3.43
3.97
p<10-3
Ageb
1.64
1.45
1.83
p<10-3
Age2b
-0.09
-0.12
-0.06
p<10-3
Litter-sizec
-1.10
-1.45
-0.76
p<lO@
Treatmentd
-0.13
-0.47
0.22
0.47
Age x Iitter-size
-0.74
-0.97
-0.51
p<10-3
Age* x litter-size
0.06
0.02
0.09
p<10-3
Age x treatment
-0.02
-0.25
0.21
0.33
Age* x treatment
0.03
-0.01
0.06
0.33
Litter-size X treatment
0.46
0.03
0.90
0.04
Age x litter-size xtreatment
0.34
0.05
0.63
0.03
Age* x litter-size x treatment
-0.05
-0.09
-0.01
0.03
a Normal approximation to the distribution of the REML estimators.
b Age and age2 were the linear and quadmtic terms for age, respectively. Age was expressed in 15day unit
and the origin was set to 15 days of age.
’ Two categories: single (reference) and multiple, the latter combining twin and triple births.
d Two categories: control (reference) and morantel.
Morantel
m-m-m
Control
o-o-o
15
30 45 60
75 90
10
Age (days)
Fig. 7. Population-averagcd live-weight growth curves açcording to litter-size and treatmcnt categories,
prçdicted from the fïxed part of the lincar-multilevcl live-weight growth mode1 of 182 Djallonke lambs bom in
10 villages in Kolda (Senegal) between October and December, 1987. The strip labels ut the top of each pane1
indicated the litter-size (single or multiple) category.
23x
R. Luncelot et ul. /Preventive Vherinury Mrdicine 46 (2000) 225-247
Table 3
Compaison of 90-day live-weights observed in different studies of Djallanke larnbs growth in West A~I%X
Keference
Country
Method
Sample size Mean
Birth type
Single
Multiple
This study (control lambs)
Senegal
LMEa
182
8.2
9.6
6.3
Agyemang et al. (1991)
T h e Gambia SMb
NA’
6.8
N A
N A
Armbruster et al. (1991)
C ô t e d’ivoire L M E
359
8.4
9.1
7.7
Charray (1986)
C ô t e d’ivoire SM
1024
9.8 NA
N A
Fall et al. ( 1983)d
Senegal
LME
360
8.1
8.7
6.5
Poivey et al. (1982)
Côte d’ivoire
OLY
293 8.6
10.0
7.1
Sumberg and Mack (1985)
Nigeria
SM
173
8.8
N A
N A
Symoens and Hxdouin (1988)
Cameroon SM
N A
11.8
13.1
9.5
Tillard (1991)
Senegal
OU
209
9.5
N A
N A
a Linear mixedeffects model.
b Sample mean.
’ Not available.
d Spline interpolation from published results.
’ Ordinxy least-squares.
Al1 the fixed effects considered important in graphical analyses were signifïcant in the
final mode1 (CI=O.~~) and those considered not important graphically (sex and its
interactions with the other explanatory variables) were not significant.
A comparison of the intercepts and slopes of the growth curves of the four important
treatment and litter-size categories (single vs. multiple litterxmorantel vs. control)
(Table 2 and Fig. 7) showed that the higher the intercept, the higher the parameter for age,
and the lower the parameter for squared-age. This latter coeffkient could be considered to
be a growth deceleration parameter.
Average daily weight gains (ADWG) in grams per day were calculated from the fixed-
effect estimates (Table 4). Initia1 differences in ADWG in the [15-30 days) age class were
associated with parameters related to linear age (agexlitter-size, agextreatment and
agexlitter-sizextreatment, see Table 2). Growth rate decreased with age. The higher the
Table 4
Averdge daily weight gains according to litter-six-by-treatment
categories, predicted by the lincar-multilevel
growth mode1 of 182 Djallonke lambs bom in 10 villages in Kolda (Senegal) between October and December,
1987 (population-averdged values)
Age ChS (dayS)
Average daiiy weight gain (grams per day)
Single lambs
Multiple lambs
Control
Morantel
Co&ol
Morantel
1 U-30)
103
103
58
78
13@-45)
91
95
53
70
[45-60)
79
86
49
63
[60-75)
67
77
45
5.5
[ 75-90)
55
68
40
48
R. Lancelot et al. /Preventive kterinary Medicine 46 (2000) 225-247
239
initial ADWG, the greater tbe decrease (i.e. growth rate decreased faster for animals with
high initial weight). These different growth decelerations were associated with
parameters related to squared-age (age x litter-size, age x treatment and age x litter-
sizex treatment, see Table 2). However, the initial treatmentxlitter-size group live-
weight rankings remained the same over a11 ages: single lambs had higher growth rates
than multiple lambs and for multiple lambs, growtb rate was higher in the morantel than
in the control group.
3.3.2. Random parameters
Random parameter estimates are listed in Table 5. The correlation was r=0.93 between
intercept weights and weights at older linear ages at the lamb level. The correlation
between the intercept and square&age weights was ~-0.60 and between the age and
squared-age weights was ~-0.26. A high positive correlation also was observed
between the intercept and linear-age weights at the herd level (~0.92).
A large overall random variation was found at the lamb-level (Table 6). The
interquartile range for the random deviation from the population mean (fixed-effect
prediction) was 2.6 kg at 90 days, and the range was 9.2 kg. At this age, population
Table 5
Random pammeters estimates for the linear-multilevel growth mode1 of 182 Djallonke lambs bom in 10 villages
in Kolda (Senegal) between October and December, 1987
Random effect
Estimate
95% confidence intervdla
Lower limit
Upper limit
Herd lrvel
S.D.(intercept)b
0.30
0.17
0.50
S.D.(age)’
0.21
0.13
0.34
Correlation (intercept, age)
0.92
0.23
0.99
S.D.(age’)”
0.04
0.02
0.05
Lamh Ievel
S.D.(intercept)
0.55
0.46
0.66
S.D.(agc)
0.33
0.27
0.40
S.D.(age2)
0.03
0.02
0.04
Correldtion(intercept,
age)
0.93
0.45
0.99
Correlation(intercept,
age*)
-0.60
-0.87
-0.06
Correlation(age, age*)
-0.26
-0.62
0.19
Order-1 nutorqressive
parameter
UI
0.6U
0.47
0.79
Hrreroscedusticiry purumeter
h
0.48
0.17
1.33
Nesidtral stnndnrd error
s
0.34
0.26
0.44
a Normal approximation to the distribution of the REML estimators.
b Standard deviation of the random effect.
’ Age was expressed in 15day units and the origin was set to 15 days of age. Age and age2 were the linear
2nd quadratic random effects for age, respectively.
240
R. Lancelot et al./Preventive Veterinary Medicine 46 (2000) 225-247
Table 6
Live-weight differences due to the rdndom effects in the linear-multilevel gmwth mode1 of 182 Djallonke lambs
born in 10 villages in Kolda (Senegal) between October and December, 1987
Age (days)
Rdndom deviation (kg)
Herd effect
Lamb effect
Overall deviation
IQRa
Rangeb
IQR
Range
IQR
Range
1 5
0.3
0.9
0.7
2.5
0.9
2.9
30
0.5
1.5
1.1
4.2
1.4
4.8
4s
0.7
2.2
1.5
5.8
1.9
6.7
60
0.9
2.9
1.8
7.0
2.3
7.8
7.5
1.2
3.7
2.1
7.8
2.8
9.7
90
1.5
4.7
2.6
9.2
3.3
11.6
a Interquartile range was the difference (due to the random effects) between the thiid and first quartiles of
live-weight devidtions from the population mean.
b Difference (due to the random effects) between the maximum and minimum values of live-weight
deviations from the population mean.
2
0
-2
2
0
-2
1.5 52 90
15 52 90
15 52 90
15 5 2 9 0
15 5290
Age (days)
Fig. 8. Deviation from population mean due to the herd-level random effects in the linear-multilevel live-weight
growth mode1 of 182 Djallonke lambs bom in 10 villages in Kokid (Senegal) between Octobcr and December,
1987. The strip label at the top of each pane1 was tbe herd identifier. Herds were ordered from the lowest (bottom
Icft panel) to the highest (top right panel) 90-day dcviation. The Upper-case letter in the top, left corner of each
pane1 stood for the treatment category: C, control; M, morantel.
R. bncelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
2 4 1
means, predicted from the fixed part of the model, ranged from 6.3 kg (multiple lambs,
control group) to 10.0 kg (single larnbs, morantel group). The herd-level variation around
the population growth curve was also large (Table 6 and Fig. 8). Treatment category was
not associated with the 90 days ranking of herd-level random deviations (Wilcoxon test:
W=648, ~28, m=17, p=O.93).
The AR(l) parameter (ai=0.66, with age expressed in 15-day periods, Table 5)
indicated that residuals’ autocorrelation persisted over at least 1 month (two 15-day
periods). The additional heteroscedasticity parameter in Table 5 (h=0.48) indicated that
the standard error of weight estimates at 90 days was about half as large as at other ages.
4. Discussion
The analysis of complex, multilevel, repeated-measures data requires that the structure
of fixed and random-effects in the mode1 be specified appropriately for correct inferences
to be made. In this example, the main study objective was to compare the growth rates of
lambs of morantel-treated vs. untreated ewes. However, a correct inference of the effect
of morantel treatment depended on: (1) identifying important fixed effects and their
interactions that could influence the morantel-treatment effect; (2) specifying the correct
random-effects structure at the multiple levels of village, herd, lamb and measurements
within lamb; and (3) specifying the correct cor-relation structure between repeated weight
measurements for individual lambs. This paper has demonstrated that graphical methods
cari play a useful role in appropriately identifying the important fixed effects and their
interactions and in specifying a logical random-effect structure to guide complex mode1
building.
In this paper, trellis graphical methods clearly and simply identified the important fixed
effects and their interactions. Al1 frxed effects highlighted by graphical analysis were
statistically significant in subsequent models and those not highlighted were insignificant.
This was the case for the lamb-sex effect and its interactions. Without graphical support
for its exclusion, we would probably have retained it as a fixed effect with interactions,
greatly increasing the rank of the fixed-effects design matrix (X) (from 12 to 24 columns).
This is important for computing time, even with simpler OLS models. With larger data
sets, calculations with 8-10 main effects and a11 their possible interactions might
overwhelm present microcomputers.
Most importantly, simple graphs highlighted that a morantel-treatment effect only
occurred in multiple-birth lambs. This lead to a mode1 formulation that identifïed a
significant morantel-treatment effect that a previous analysis, using cross-sectional
analysis-of-variante (ANOVA), had missed (Tillard, 1991). At the time of this first
analysis, no software package was widely available to fit multilevel models. Because
of the lack of power of the cross-sectional ANOVA approach, litter-size-by-treatment,
and age-by-litter-size-by-treatment interactions were not found to be significant
CJDO.05).
Graphical methods also proved very useful at identifying random-effects structures or
potential problems with their estimation at multiple levels. This is crucial when analysts
are confronted to new kinds of data without preliminary information on correlations
242
R. Lancelot ef al. /Preventive krerinary Medicine 46 (2000) 225-247
between random effects. From our experience, it is not straightforward to choose an
appropriate (explanatory and parsimonious) structure for the random-effects variance-
covariance matrix. Cubic (or higher-order) polynomials might be necessary to mode1
live-weight growth over a longer period of time (e.g. from birth to 1 year). The number
of variance-covariance parameters in an unstructured matrix, 2, is large (10 for each
hierarchical level with a cubic polynomial) and convergence problems are likely to
occut In such situations, graphical analysis might be the only reasonable (or even
possible) way to specify a plausible and numerically stable random structure. In this
analysis, initial graphical assessments were consistently confïrmed in subsequent models.
At the herd level, growth pattems varied and this Was reflected in significant herd
random intercepts and slopes. At the lamb level, trellis plots revealed more-complex
correlation patterns (leading to model-building strategies starting with a complex
correlation structure). For individual weight measurements within lamb, graphical
methods to identify important correlations at different lags and differences in residual
variation at different ages also provided useful information. In our opinion, the graphical
analyses used were a powerful tool in guiding mode1 building and in effectively
displaying important results.
Dam might have been considered as a grouping level in this analysis. As a matter of
fact, lambs bom from multiple litters might share common and unobserved genetic,
nutritional and health features, inducing cor-relation among measurements from the same
litter. This level was tested (results not shown here) but no significant random effect nor
covariance parameter could be identifïed. This is probably because most of births were
single or of small size (182 lambs from 147 ewes and 147 litters: 104 single, 68 twin and
10 triple lambs) and/or the litter-size fixed effect accounted for within-litter correlation.
With more-prolific species such as West-African-Dwarf goats (for which we have similar
data sets) over a longer observation period (with several litters per dam), a dam effect
might be more important.
The 90-day live-weight gain results in this study were within the range of observations
made in other studies in the region (Table 3). The lowest weight gains were observed in
Gambia (Agyemang et al., 1991); however, few details were given of the study methods.
The largest weight gains were observed in northwest Cameroon (Symoens and Hardouin,
1988) in an experimental herd under improved management (housing, supplementary
feeding, herd health). Thus, their results were not likely to be representative of the usual
growth performance in the local farming system. The other estimates in Table 3 were
similar to those reported in this study - particularly the magnitude of the multiple birth
(litter-size) effect.
The effect of litter-size on weight gain was large. Multiple larnbs were 1.1 kg lighter
than single lambs at 15 days of age and their growth rate was lower (regardless of age and
treatment group). Low growth rate for multiple lambs indicated that their nutritional
needs were not met by the ewes’ milk yields. At the end of the study period, single lambs
bom from morantel-treated ewes had a slightly greater growth rate than lambs bom from
control ewes (Fig. 7). However, the important effect of morantel treatment was in
multiple-birth lambs. At 90-day of age, multiple lambs bom from treated ewes were
1.3 kg heavier, on average, than those bom from control ewes (7.6 vs. 6.3 kg,
respectively, i.e. 22% more). The morantel effect was thought to be due to a better
R. Luncelot et al:/Preventive kterinary Medicine 46 (2000) 225-247
243
milk yield for the dewormed ewes - crucial for ewes with multiple lambs (Theriez,
1984). This study took place just after the rainy season, when nutritional and health
conditions were at their best with abundant forage resources and trop residues and low
parasitic pressure. Morantel treatment might have an even greater effect in lambs born
during the second lambing peak from March to May, when ewe nutritional stress would
be greatest. Because multiple births were frequent in Kolda (43% of the lambings in this
sample), the use of morantel for deworming ewes during the rainy season should have an
important effect on flock productivity. A thorough assessment of the benefit of
deworming should take into account its effect on other. production parameters (e.g.
mortality, fertility, prolificacy) as well as the cost of treatment. This was done in a
companion study (Lesnoff et al., 2000). Using a seasonal population-dynamics modeling
approach (steady-state model), a positive effect of deworming was found, with the
financial benefit-cost ratio estimated at 3.7 11.9, 5.41 (95% confidence interval in
brackets).
However, the persistency of this positive morantel-effect on growth at older ages needs
to be assessed in further fïeld studies. Published data are available with other drugs -but
to our knowledge, not with morantel for Djallonke sheep.
At the herd and lamb levels, intercept and linear-age random effects were highly
correlated. The persistence of this cor-relation at older ages should be checked to
determine whether early weights are as good predictors for mature weight as for 90-day
weights. The major consequence would be for the selection of females. In an experiment
with Scottish Blackface sheep, Gunn (1977) demonstrated that ewe mature size, sexual
precocity and lifetime productivity cari be affected by insuffrcient early growth. However,
in Nigeria, an analysis of station records of Yankasa ewes showed no association between
sexual precocity and litter-size, birth weight and weaning weight (Osuhor et al., 1997).
However, these are not strong data (observational rather than experimental and from a
specific station).
The identification of herd and lamb random effects provide important targets for
specific studies to improve lamb growth. The presence of herd random effects is
potentially related to differences in management, feeding and herd-health practices.
T~US, studies contrasting these features in flocks with low and high growth rates
should be considered. Lamb random effects might be associated with genetic
differences, ewe nutritional and health status, lamb diseases and other factors at the
lamb level.
Graphical methods are useful in highlighting the best- and worst-performing flocks or
herds. A quick examination of Fig. 8 shows that 20% of herds (nine out of 45, bottom
line) were responsible for rnost low weights and approximately 20% for most high
weights (top line). Statistical methods that ignore these herd differences are likely to
rnake incorrect inferences on fixed effects - particularly with respect to too small
variante estimates (McDermott et ai., 1994). There are particular benefrts in better
specifying the fixed and random-effects structures with graphical methods. Random-
effect deviations from population averages are intuitively appealing in such multilevel
field studies - and if correctly specified, are more powerful at detecting differences than
rnarginal (e-g. Zeger and Liang, 1986) or generalized-least-squares (e-g. Box et al., 1994)
models (van der Leeden et al., 1996).
744
R. /~~nc&t er al. /Prevertlive Veteritm-y Medicitw 46 (2000) 225.247
Another important modeling issue to which graphical tools cari be applied is in the
decision of how parsimonious or complex a mode1 should be. For example, in these data,
autocorrelation pattems were complex. Even an unstructured L$ matrix (one parameter
for each lag pair, for a total of 21 correlation parameters) could not remove residual
autocorrelation. Graphical assessment (Fig. 2) showed that a small proportion of lambs
(less than 10%) grew irregularly and poorly and that their growth was not easily modeled
by age and age-squared random effects. Deleting these poor-growing lambs to improve
mode1 fit is not an option. Heroic measures to better mode1 their growth is probably not
wise either. The use of a simple AR(l) mode1 to parsimoniously explain residual
autocorrelation for most of the lambs, the recognition that poor-growing lambs bave hard-
to-mode1 growth curves (for which specific studies could be initiated if considered
important) and the use of graphical methods to investigate when poor mode1 fit may
influence inferences about fixed and random effects seems a sensible compromise.
Structured graphical analysis with trellis graphs and other recently developed methods are
important tools in efforts to effectively summarize and test biological processes with
sufficient and not excessive complexity.
5. Conclusion
The graphical methods presented here are an important adjunct to the consideration of
linear-multilevel growth models. Users cari specify relevant starting points for the fixed
and the random effects, check mode1 goodness of fit and assess distributional assumptions
within the same computing and graphie environment. Until recently, this process was
cumbersome and required different software packages to perform different steps
(Mazumdar et al., 1999). The use of graphical methods in mode1 building and assessment
should greatly enhance the ability of veterinary epidemiologists to explore different
random structures and to effectively communicate their fïndings.
In this example, graphical analyses greatly aided in highlighting important fixed and
random effects. Of particular importance was the positive effect of ewe deworming on the
pre-weaning live-weight growth of multiple-birth lambs. In addition, the relationships
between early weight measures and subsequent growth and the large variability of growth
between herds were nicely demonstrated by graphical anal,ysis of correlation patterns and
herd random effects. In future, we foresee that graphical tools Will become more user-
friendly, allowing for interactive outputs to explore specifrc questions of both researchers
and their clients. This should allow for more their broader application and the better
integration of results from individual studies into a more general livestock systems
context.
Acknowledgements
This work was funded by a joint ISRA-CIRAD research programme f’uthofogie er
p-oductivite’ des petits ruminants dans les systèmes traditionnels (Small-ruminants
production and health in traditional farming systems). We thank Olivier and Brigitte
R. -celot et af. /PreVentiVe Vetetinary Medicine 46 (2000) X5-247
245
Faugère for providing us with data, and ISRA, Senegalese Veterinary Services and
farmers for their co-Op=tion. We also thank José Pinheiro (Lucent Technologies,
Murray HiIl, NJ) and Douglas Bates (University of Wisconsin, Madison, WI) for their
help in handling and understanding linear-multilevel models, and for providing us wit,h an
early version of their book. The two anonymous reviewers provided helpful comments on
the fiist draft of this paper.
References
Abassa, K.I?, Pessinaba, J., Adeshola-Ishola,
A., 1992. Croissance pi-sevrage des agneaux Djallonké au centre
de Kolokopé (Togo). Rev. Élev. Méd. Vét. Pays Trop. 45 (l), 49-54.
Agyemdng, K., Rawlings, P., Clifford, D., Bojang, N., Samba, A., 1991. Productivity and health parameters of
mal1 ruminants in villages of the Gambia. Bull. Anim. Hlth. Prod. Afr. 39, 129-135.
Akaike, H., 1974. A new look at statistical mode1 identification. IEEE Trans. Automatic Control AU-19, 716-
722.
-bruster, T., Peters, K.J., Metz, T., 1991. Sheep production in the humid zone of West Africa. II. Growth
performance and liveweight of sheep in improved and h-aditional production systems in Côte d’ivoire. J.
Anim. Breed. Genet. 108, 210-219.
Becker, R.A., Cleveland, W.S., Shyu, M.-J., 1996. The visual design und control of trellis display. J. Comput.
hphicti Statist. 5 (2), 123-155.
Box, G.E.P., Jenkins, G.M., Reinsel, G.C., 1994. Time Series Analysis: Forecasting and Contml, 3rd Edition.
Prentice-Hall, Englewood Cliffs, NJ, 598 pp.
Charray, J., 1986. Performances des brebis Naines de I’Afrique de I’Ouest entretenues suivant deux rythmes
différents d’accélération de h reproduction, Rev. Élev. Méd. Vét. Pays. Trop., 39 (1) 15 l-160.
Clément, V., Poivey, J.P., Faugère, O., Tillard, E., Lancelot, R., Gueye, A., Richard, D., Bibé, B., 1997. Etude de
la variabilité des caractères de reproduction chez les petits ruminants en milieu d’élevage traditionnel au
Sénégal. Rev. Élev. Méd. Vét. Pays Trop. 50 (3), 235-249.
Cleveland, W., 1993. Visualizing Data. Summit. Hobart Press, 360 pp.
Cleveland, W., 1994. The Elements of Graphing Data. Summit. Hobart Press, 297 pp.
Davidian, M.M., Giltinan, DM., 1995. Nonlinear Models for Repeated Measurement Data. Chapman & Hall,
London, 359 pp.
Diggle, P.J., 1990. Tiie Series. A Biostatistical Introduction. Clarendon Press, New York, 257 pp.
Diggle, P.J., Liang, K.-Y., Zeger, S.L., 1994. Analysis of Longitudinal Data. Clarendon Press, Oxford,
253 pp.
Fall, A., Diop, M., Sandford, J., Gueye, E., Wissocq, Y.J., Durkin, J., Trail, J.C.M., 1983. Étude sur la
productivité moutons Djallonké au Centre de Recherches Zootechniques de Kolda, au Sénégal. 2. Poids
corporel, productivité des brebis et du troupeau. Rev. Élev. Méd. Vét. Pays Trop. 36 (3), 283-289.
Paugère, O., Faugère, B., 1986. Suivi de troupeaux et contrôle des performances individuelles des petits
ruminants en milieu traditionnel africain. Aspects méthodologiques. Rev. Élev. Méd. Vét. Pays Trop. 39 (l),
29-40.
Wugère, O., Faugère, B., 1993. Panurge: suivi individuel dans les systèmes d’élevage traditionnel. Maisons-
Alfort (Frdnce), CIRAD-EMVT, 339 pp.
Faugère, O., Dockès, A.-C., Perret, C., Faugère, B., 1990. L’élevage traditionnel des petits ruminants au Sénégal.
1. Pratiques de conduite et d’exploitation des animaux chez les éleveurs de la région de Kolda. Rev. Élev.
Méd. Vét. Pays Trop. 43 (2), 249-259.
Faugère, O., Tillard, E., Faugère, B., 1990. Prophylaxie chez les petits ruminants au Sénégal: régionalisation
d’une politique nationale de protection sanitaire. In: Rey, B., Lebbie, S., Reynolds, L. (Eds.), First Biennial
Conference of the Africun SKndll Ruminant Research Network, Nairobi, INCA, pp. 307-314.
Fri&che, T., Kaufmann, J., Mister, K., 1993. Parasite spectrum and seasonal epidemiology of gastrointestinal
nematodes of small ruminants in the Gambia. Vet. Parasitol. 49, 271-283.
R. hcelrzr et al./Preventive Veterinaty Medicine 46 (2000) 225-247
Gol&t&,, H., 1986. Multihd mixed linCX mode1 andySiS USing iterdtive generdlized least squares. Biomet&+
73 (1). 43-56.
Goldst&, H., 1995. MultiIevel Statistical Models. Arnold, London, 178 pp.
cEen, L-E., ~~~~~ E., Morgan, K.L., 1998. A multi-level mode1 of data with repeated measures of the effect
of bb d&&fja on weight. Prev. Vet. Med. 36 (2), 85-94.
Gc~, Y-T., McDermott, J.J., Schukken, Y.H., Hertl, J.A., Eicker, S.W., 1999. Analysis of con-ekited continuous
reps observations: modelling the effect of ketosis on milk yield in dairy cows. Prev. Vet. Med. 39, I37-
153.
Gunn, R.C., 1977. The effects of two nutritional envimnments from 6 weeks pre partum to 12 months of age on
lifetime performance and reproductive potential of Scottish Bhu&face ewes in two adult environments.
Anim. Prod. 25, 1.55-164.
Kraft, I., de L.eeuw, J., van der Leeden, R., 1994. Review of fïve rnultilevel analysis pmgrams : BMDP-SV,
GENMOD, HLM, h4L3, VARCL. The Am. Statist. 48 (4), 324-335.
Laird, N.M., Ware, J.H., 1982. Random-effectb: models for longitudinal data. Biometrics 38, 963-974.
Lancelot, R., Faye, B., Ju&s, X., Ndiaye, M., Pérochon, L., Tillard, E., 1998. La base de données BAOBAB: un
outil pour modéliser la production et la santé des petits ruminants dans les systèmes d’élevage traditionnels
au Sénégal. Rev. Élev. Méd. Vét. Pays Trop. 51 (2), 135-146.
Lee, JC., 1991. Tests and mode1 selection for the general growth curve model. Biometrics 47, 147-159.
Lesnoff, M., Lancelot, R., Tillard, E., Dohoo, LR., 2000. A steady-state appmach of benefit-cost analysis with a
periodic Leslie-m&% model: presentation and application to the evduation of a sheepdiseases preventive
scheme in Kolda Senegal. Prev. Vet. Med. 46, 113-128.
Littel, R.C., MilIiken, C.A., Stroup. W.W., Wolfimger, R.D., 1996. SAS System for Mixed Models. Cary, NC,
663 PP.
Mazumdar, S., Begley, A.E., Houck, P.R., Yang III, Y., Kupfer, D.J., 1999. Residual analysis in random
regrcssion using SAS and S-PLUS. Comput. Meth. Prograns Biomed. 58, 281-282.
McDermott, JJ., Schukken, Y.-H., Shoukri, M.M., 1994. Study design and analytic methods for data collected
from clusters of animals. Prev. Vet. Med. 18, 17.5-191.
McDermott, J.J., O’Callaghan, C., Kadohira, M., 1997. Experiences in multi-level models for assessing the
distribution of animal diseases. Épidémiol. Santé Anim. 31/32, 13.23.1-13.23.3.
Orji, B.I., Steinbach, J., 1981. Post-weaning growth and development of Nigerian dwarf sheep. Trop. Anim.
Hlth. Pmd. 13, 101-106.
Osuhor, C.U., Osinowo, O.A., Nwagu, B.I., Dennar, F.O., Abdullahi-Adee,
A., 1997. Factors affecting age at
fmt lambing in Yankasa ewes. Anim. Reprod. Sci. 47, 205-209.
Patterson, H.D., Thompson, R., 1971. Recovery of inter-block information when block sizes are unequaI.
Biometrika 58, 545-554.
Pinheiro, J.C., Bates, D.M., 1996. Unrestricted parameterizations
for variance~ovariance
matrices. Statist.
Comput. 6, 289-296.
Pinheiro, J.C., Bates, D.M., 2000. Mixed-effects Models in S and S-Plus. Springer, New York, 598 pp.
Poivey, J.P., Landais, E., Berger, Y., 1982. Etude et amélioration génétique de la croissance des agneaux
Djallonké. Résultats obtenus au Centre de Recherches Zootechniques de Bouaké (Côte d’ivoire). Rev. Élev.
Méd. Vét. Pays Trop. 35 (4), 421-433.
Puyalto, C., Sanaa, M., N’Djoya, A., Planchenault, D., 1997. Factors influencing daily weight gain on Young
zebu in Garoua district, Camemon. Épidémiol. Santé Anim. 31/32, 02.03.1-02.03.3.
Rao, C.R., 1965. The theory of least squares when the parameters are stochastic and its application to t.he
analysis of growth curve. Biometrika 52 (3/4), 447-458.
Sumberg, J.E., Mack, S.D., 1985. Village production of West African Dwxf go& ad sheep in Nigeria
Trop. Anim. Hlth. Prod. 17, 135-140.
Symoens, C., Hardouin, J., 1988. Le mouton Djallonké en élevage extensif dans le Nord-Ouest Cameroun. Rev.
Élev. Méd. Vét. Pays Trop. 41 (4), 449-458.
van der Leeden, R., Vrijburg, K., de Leeuw, J., 1996. A review of two different approaches for tbe andysis of
growth data using longitudinal mixed linear models: comparing hierarchical linear regression (ML3, HLM)
and repeated maures designs with stmctured covtiance matrices (BMDPSV), The Statist. Software
Newslett. 20, 583-605.
R. Lmcelot et al. /Preventive Veterinary Medicine 46 (2000) 225-247
247
Venables, W.N., Ripley, B.D., 1999. Modem Applied Statistics with S-Plus. Springer, New York, SO1 pp.
Theriez, M., 1984. Influence de l’alimentation sur les performances de reproduction des ovins. In: Neuvièmes
Journées de la Recherche Ovine et Caprine, Paris, INRA-ITOVIC, pp. 294-326.
Tillard, E., 1991. Évaluation technico-économique des effets des prophylaxies en milieu villageois chez les
petits ruminants de la région de Kolda (Sénégal). Thèse Méd. Vét. ENVA. Faculté de Médecine de Créteil,
172 pp.
Tuah, A.K., Baah, J., 1985. Repmductive performance, pre-weaning growth rate and pre-weaning lamb
mortality of Djallonke sheep in Ghana. Trop. Anim. Hlth. Prod. 17, 107-113.
Zeger, S.L., Liang, K.Y., 1986. Longitudinal data analysis for discrete and continuous outcomes. Biometrics 42,
121-130.