A sinusoidal function and the Nelder-Mead simplex algorithm applied to growth data from broiler chickens

A sinusoidal function and the Nelder-Mead simplex algorithm applied to growth data from broiler chickens H. Darmani Kuhi,∗,1 A. Shabanpour,∗ A. Mohit,∗ S. Falahi,† and J. France‡ ∗

optimization method based on a direct search method instead of using gradient-based techniques. The results of this study show that both the sinusoidal function and the direct search method precisely describe the growth dynamics of broiler chickens. Fitting the growth functions to different data profiles nearly always led to the same or less maximized log-likelihood values for the sinusoidal equation, which is an indication of its superiority in describing growth data from broiler chickens.

ABSTRACT There has been much recent interest in mathematical developments for the analysis of growth in poultry. In this paper, we present a sinusoidal function to describe the evolution of growth as a function of time based on real life experiments. The function was evaluated with regard to its ability to describe the relationship between body weight and age in broilers and was compared to 4 standard growth functions: Gompertz, logistic, Lopez, and Richards. In order to estimate the model parameters, we adopted a global

Key words: Growth functions, sinusoidal equation, Nelder-Mead algorithm, broilers 2018 Poultry Science 97:227–235 http://dx.doi.org/10.3382/ps/pex299

INTRODUCTION

should describe data well and allow biologically and physically meaningful traits to be determined (France et al., 1996). Today, different methods are available to the analyst for optimizing nonlinear functions in order to fit various mathematical prediction models to experimental data. These are termed gradient-based algorithms. They are necessarily iterative in nature and need to be provided with initial values for the model parameters. These values are usually chosen either by making an educated guess based on the structure of the model, or by using parameters estimated from previous studies. Many of these methods (e.g., Trust Region, Newton-Raphson with line search, quasi-Newton, double-dogleg, and conjugate gradient methods) use first-order or sometimes second-order derivatives to determine a search direction to improve the value of the objective function (Fletcher, 1981; Sorensen, 1982; Dennis and Schnabel, 1983; Luenberger, 1984; Nocedal and Wright, 1999). However, gradient-based algorithms (e.g., Gauss-Newton, Levenberg-Marquardt) are totally deterministic, and when applied to problems affected by noise, whether it be error in measurement or uncertainty in prediction, they are either unable to reach an optimum or they may reach a local optimum (Young, 1976, pp. 517–518). An incorrect choice of starting values can make the algorithm head towards a local optimum, a place where many sets of parameters produce a similar likelihood. In these cases, for example in all Newtonian methods, the problem becomes highly

Poultry industries face various decisions in the production cycle that include nutrient and mineral supply to birds, cost and type of feed, and a range of bird health, welfare, and environmental issues that affect the profitability of an operation. Prediction of growth and identifying the times of maximum growth rate and when the birds are ready for sale are important factors that contribute to the profitability of poultry operations. Traditionally, mathematical equations called growth functions have been used to relate body weight (BW) to age of the bird or cumulative dietary feed intake (McCance, 1960; Lister et al., 1966; Fitzhugh, 1976; Darmani Kuhi et al., 2003, 2004; Faridi et al., 2014, 2015). They can also be used to determine the efficiency of nutrient utilization, which is the derivative of the relationship between BW and dietary nutrient intake and as response functions to predict daily energy, protein, and amino acids requirements for maintenance and growth (France et al., 1989; Darmani Kuhi et al., 2009, 2011, 2012). Such applications have been reported in animal science for nearly a century (for review see Dumas et al., 2008; Darmani Kuhi et al., 2010). A useful growth function

 C 2017 Poultry Science Association Inc. Received September 26, 2016. Accepted September 15, 2017. 1 Corresponding author: darmani [email protected]

227

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

Department of Animal Science, Faculty of Agricultural Sciences, University of Guilan, Rasht, Iran; † Department of Mathematics, Salman Farsi University of Kazerun, Kazerun, Iran; and ‡ Centre for Nutrition Modelling, Department of Animal Biosciences, University of Guelph, Guelph ON, N1G 2W1, Canada

228

DARMANI KUHI ET AL.

MATERIALS AND METHODS Sinusoidal Model

W = a1 + a2 sin



−

n  π t − a3 2k

(1)

where W is live weight (g) at age t (d) and a1 , a2 , a3 , k and n are parameters. The above sinusoidal equation provides a flexible growth function which can describe sigmoidal and diminishing returns behavior. It does well at describing sigmoidal patterns (Figure 1) exhibiting

faster early growth and a fairly low but variable point of inflection. This function also can describe a wide range of hyperbolic shapes when there is no point of inflection. Since the sinusoidal function is periodic, we only need consider its behavior over a particular interval, e.g., the green line in Figure 1, where [1] a1 : the vertical offset (height of the baseline) is the average value of the growth coefficient during an oscillation. When a1 is changed, the basic sinusoidal function is shifted vertically by a1 units. [2] a2 : the amplitude (the height of each peak above the baseline) of curve is the overall rate parameter of the entire growth process. Multiplication of the function by 2 doubles the amplitude, or distance it oscillates up and down. While the original curve oscillated between −1 and +1 based on the baseline, the new curve oscillates between −2 and +2. In general, a2 sin(x) has amplitude a2 based on the baseline. Note that |a1 | + |a2 | is an approximation to the final weight (weight at maturity). [3] a3 : is the phase shift (the horizontal offset of the base point). By changing a3 , the amplitude and period of the function appears to be the same, but the graph of the function has been shifted horizontally by a3 units. This horizontal shift is generally called the phase shift. The phase shift a3 tells us when the curve first crosses the time axis as it ascends. [4] L: the period of the sinusoidal function, which is twice the length of growing period. [5] ω : the angular frequency of the function, given by ω = 2Lπ . As we can see, L = −4k and so k is −L/4. Parameters a1 , a2 , a3 , and k are all real numbers while n is real but positive. Point of inflection determines the time at which the absolute growth rate is maximum and it plays an

Figure 1. Graph of the sinusoidal function showing its fit to the data of Hancock et al. (1995). Where a1 is the height of the baseline; a2 is the height of each peak above the baseline; a3 is the horizontal offset of the base point, where the curve crosses the baseline as it ascends initial weight and L is the period of the sinusoidal function.

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

identified and then inversion of the Hessian matrix fails (Juretig et al., 2012). In addition to these gradient-based algorithms, there is a class of techniques known as direct search methods. In contrast to other optimization techniques which require derivatives of the objective function to determine a search direction, a direct search method relies only on the value of the objective function at a set of points (Hooke and Jeeves, 1961; Torczon, 1993). A direct search method does not require any derivative information, which makes it suitable for problems with non-smooth functions. With a direct search method, the calculations are simple and few adjustable parameters need to be supplied. This class of methods is widely used to solve for parameter estimates and similar statistical problems where the function values are uncertain or subject to noise. It can also be used for problems involving discontinuous functions, which occur frequently in statistics and experimental mathematics. The objective of the present study is to introduce a sinusoidal function into poultry science by applying it to temporal growth data from broiler chickens using a specific direct-search method (e.g., the Nelder-Mead algorithm) for parameter estimation.

229

SINUSOIDAL FUNCTION FOR GROWTH IN BROILERS Table 1. Some properties of the different reference growth functions considered. Time at inflection point (t∗ )

Growth function

Functional form1

Gompertz

W = W 0 exp[(ln

logistic

W =

W0Wf W 0 +(W f −W 0 )e−b t

b

Lopez

W =

(W 0 K n +W f t n ) (K n +t n )

1/ n K [ nn −1 +1 ]

Richards

W =

1

Wf W0

)(1 − e−b t )]

W0Wf

1/ n (W 0 n +(W f n −W 0 n ) exp(−b t ))

W

b

1

[ln(ln( Wf0 )]

1

ln(

1

b

ln(

0.368W f

W f −W 0 ) W0

W n −W 0n f

n W 0n

Weight at inflection point (W∗ )

0.5W f [(1+ n1 )W 0 +(1− n1 )W f ] 2

Wf

)

(n +1)1/ n

important role in growth curve studies across disciplines. The sinusoidal equation has a flexible inflection point (W∗ , t∗ ) that occurs at:

f (x1 ) ≤ f (x2 ) ≤ ... ≤ f (xn +1 )

2ka3 W ∗ = a1 ; t∗ = − π provided n > 2. This can be seen by equating zero. Differentiating Equation (1) twice gives

[1] Order. Order according to the values at the vertices:

d2 W dt 2

to

 π 2 d2 W = na θn −2 [(n − 1) cos (θn ) −nθn sin (θn )] ; 2 dt2 2k π θ = − t − a3 (2) 2k

Equation (2) equates to zero when θ equals zero, i.e., t = −2ka3 /π , provided n > 2. To examine the efficacy of the sinusoidal model in describing the growth pattern in broiler chicks, 4 classical growth functions, e.g., Gompertz, logistic, Lopez, and Richards (Thornley and France, 2007, Chapter 5), were used as reference in this study (Table 1).

Nelder-Mead Algorithm The Nelder-Mead algorithm is simplex-based direct search method (Nelder and Mead, 1965). A simplex S in n-dimensional real space (Rn ) is defined as the convex hull of n + 1 vertices x1 , ..., xn +1 ∈ Rn , where Rn is the set of all n dimensional vectors. The method begins with a set of n + 1 points x1 , ..., xn +1 ∈ Rn which are considered as the vertices of a working simplex S, and the corresponding set of real function values at the vertices, e.g., f (xj ) = 0; j = 1, ..., n + 1. The method then performs a sequence of transformations of S, aimed at decreasing the function values at its vertices. At each step, the transformation is determined by computing one or more test points, together with their function values, and by comparison of these functions values with those at the vertices. This process is terminated when S becomes sufficiently small in some sense, or when the function values f(xj ) are close enough. An outline of the steps involved follows.

[2] Centroid. Calculate x0 , the centroid of all points except xn+1 :

x0 =

n 



xj

n

j =1

[3] Reflection. Compute the reflected point:

xr = x0 + α (x0 − xn +1 ) If the reflected point is better than the second worst but not better than the best, i.e., f (x1 ) ≤ f (xr ) ≤ f (xn ) then obtain a new simplex by replacing the worst point xn+1 with the reflected point xr and go to step 1. [4] Expansion. If the reflected point is the best point so far, f (xr ) ≤ f (x1 ) then compute the expanded point:

xe = x0 + γ (x0 − xn +1 ) If the expanded point is better than the reflected point, f (xe ) ≤ f (xr ) then obtain a new simplex by replacing the worst point xn+1 with the expanded point xe , and go to step 1, else obtain a new simplex by replacing the worst point xn+1 with the reflected point xr and go to step 5. [5] Contraction. Here, it is certain that f (xr ) ≥ f (xn ). Compute contracted point:

xc = x0 + ρ (x0 − xn +1 ) If the contracted point is better than the worst point, i.e., f (xc ) ≤ f (xn +1 ), then obtain a new simplex by replacing the worst point xn+1 with the

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

W is body weight; t is time; Wf is final weight, W0 is initial weight, and b, K and n are parameters.

230

DARMANI KUHI ET AL. Table 2. Data sources1 used in the study. Growth phase d

M2

F2

Strain

Considerations

Leeson and Summers (1980) Plavink and Hurwitz (1983) Hancock et al. (1995)

7-77 0-70 1-48

6 1 6

6 1 6

Carcass analysis study Carcass analysis study Genetic selection study

Govaerts et al. (2000)

7-70

1

1

Broilers White rock Hubbard, Hybro, Ross 608, Ross 688, Ross 708 and Ross 786 Slow and fast growing broilers

Source

1 2

Quantitative feed restriction study

The growth data are the average values presented in the referenced papers. Male = male, F = female. The numbers under M and F refer to the number of data profiles from a specific data source.

Figure 3. Plots of body weight (g) against age (d) showing the fit of different growth functions to female data profile of Hancock et al. (1995) (Ross 788).

contracted point xc and go to step 1. Else go to step 6. [6] Reduction. For all but the best point, replace the point with:

xi = x1 + σ (xi − x1 )

for all i ∈ {2, 3, ..., n + 1} . Go to step 1. Note that α, γ, ρ and σ are real scalars, and are the reflection, expansion, contraction, and shrink coefficients, respectively. Standard values are α = 1, γ = 2, ρ = −0.5 and σ = 0.5. Model parameters were estimated using the MATLAB nlmefit (MATLAB 8.0, The MathWorks Inc., Natick, MA, 2012) function that

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

Figure 2. Plots of body weight (g) against age (d) showing the fit of different growth functions to male data profile of Hancock et al. (1995) (Ross 708).

231

SINUSOIDAL FUNCTION FOR GROWTH IN BROILERS

Table 3. Estimated growth parameters for given male data profiles using the different growth functions. Standard error estimates for parameters are given in brackets. Growth parameters1 Function

W0

Wf

b

K

n

a1

a2

a3

Statistical Criterion2 RMSE2

R2

∗

Leeson and Summers, 1980 Gompertz Logistic

Richards Sinusoidal

4669.5 (23.3) 3630.1 (152.9) 5417.54 (1263.8) 3786.3 414.7) 3635.3

0.036 (0.004) 0.077 (0.005)

52.5 (19.8) 128.6 (29.0) 39.7 (94.9) 15.3 (46.2) 45.3

4368.1 (435.7) 3364.3 (226.7) 6829.1 (2937.5) 5653.3 (2679.5) 3953.7

0.034 (0.004) 0.072 (0.007)

57.3 (5.31) 160.9 (18.8) 90.1 (26.8) 66.1 (12.4) 53.8

5795.2 (97.9) 4672.0 (122.1) 7234.2 (459.3) 5645.4 (199.9) 4617.2

0.036 (0.001) 0.073 0.003)

51.6 (4.56) 99.7 (13.1) 82.2 (21.8) 45.4 (12.9) 49.6

4190.7 (179.2) 2751.3 (151.5) 6486.1 (1357.8) 4578.4 (879.7) 3058.2

0.039 (0.001) 0.094 0.006)

58.7 (12.1) 0.068 (0.023)

−52.7 (0.005)

2.31 (0.399) 0.706 (0.560) 0.826 (0.051)

1645.5 (258.3)

1990.0 (99.5)

95.5 (0.237)

80.49

99.50

73.09

99.59

87.94

99.41

74.06

99.58

72.72

99.60

61.87

99.60

84.00

99.26

60.54

99.62

60.62

99.62

62.52

99.59

22.03

99.98

68.66

99.79

32.48

99.95

22.17

99.98

12.07

99.99

9.76

99.98

31.81

99.82

7.87

99.99

8.55

99.99

8.50

99.99

Plavink and Hurwitz, 1983 Gompertz Logistic Lopez Richards Sinusoidal

81.3 (32.9) 0.021 (0.016)

−51.8 (0.005)

1.90 (0.386) −0.341 (0.396) 1.10 (0.041)

1998.8 (328.3)

−1954.9 (283.0)

92.6 (0.135)

Hancock et al., 1995 (Ross 788) Gompertz Logistic Lopez Richards Sinusoidal

64.7 (3.78) 0.038 0.003)

−48.1 (0.0003)

2.19 (0.100) 0.067 (0.085) 0.969 (0.176)

2326.2 (27.9)

2291.0 (23.5)

95.7 (0.046)

Govaerts et al., 2000 (Fast-H) Gompertz Logistic Lopez Richards Sinusoidal 1 2

69.1 (12.1) 0.034 (0.010)

−40.3 (0.003)

2.01 (0.139) −0.090 (0.168) 0.967 (1.79)

1546.3 (133.6)

1511.9 (100.2)

64.3 (0.34)

W0 is initial weight, Wf is final weight, and a1 , a2, a3, b, K, and n are parameters. RMSE = Root mean square error, and R2 = Adjusted r squared.

allows fitting a nonlinear mixed-effects model. By default, nlmefit fits a model in which each parameter is the sum of a fixed and a random effect, and the random effects are uncorrelated (their covariance matrix is diagonal). The maximized log-likelihood (MLE), estimated parameters of the error variance (Errorparam), Akaike information criterion (AIC), Bayesian information criterion (BIC) and estimated error variance (MSE) were used to evaluate the general goodness-of-fit of each model to the different data profiles. The higher values of MLE are preferred. AIC and BIC are criteria for model selection among a finite set of models; the model with the lowest AIC and BIC is preferred. These

are based, in part, on the likelihood function. When fitting models, it is possible to increase the likelihood by adding parameters, but doing so may result in overfitting. Both BIC and AIC attempt to resolve this problem by introducing a penalty term for the number of parameters in the model; the penalty term is larger in BIC than in AIC.

Data Source and Statistical Analysis Twenty-eight time course profiles from 4 different studies (Table 2) were used in this study to investigate the relationships between BW and age in different male

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

Lopez

54.7 (454.9) 127.6 (22.1) 119.792 (74.1) 110.3 (43.1) 36.7

232

DARMANI KUHI ET AL.

Table 4. Estimated growth parameters for given female data profiles using the different growth functions. Standard error estimates for parameters are given in brackets. Growth parameters Function

W0

Wf

b

K

n

a1

a2

a3

Statistical Criterion2 RMSE

R2

46.71

99.73

70.09

99.39

45.41

99.75

44.56

99.76

41.23

99.79

44.38

99.69

62.37

99.39

44.76

99.69

45.31

99.68

38.84

99.76

23.77

99.96

76.22

99.59

15.88

99.98

15.86

99.98

24.59

99.96

9.99

99.97

25.49

99.81

10.15

99.97

9.86

99.97

9.86

99.97

Leeson and Summers, 1980 Gompertz Logistic

Richards Sinusoidal

4215.1 (377.1) 3102.3 (195.3) 9286.1 (5180.2) 5883.2 (2717.7) 3725.75

0.031 (0.003) 0.068 (0.006)

36.7 (13.0) 101.6 (21.1) 65.8 (62.8) 25.9 (35.0) 76.6

3128.5 (198.4) 2570.5 (120.8) 3877.9 (717.1) 3296.3 (604.8) 3566.1

0.040 (0.003) 0.078 (0.007)

71.2 (7.10) 162.9 (24.0) 65.0 (14.1) 43.2 (11.1) 40.1

4561.1 (96.7) 3726.3 (134.1) 6419.9 (310.8) 5050.0 (235.5) 3699.2

0.035 (0.001) 0.070 (0.004)

52.3 (5.43) 96.5 (12.9) 76.4 (24.7) 43.7 (15.3) 36.43

3135.7 (146.6) 2205.1 (119.0) 4829.8 (1102.4) 3505.0 (740.5) 2442.313

0.042 (0.002) 0.093 (0.006)

122.9 (67.0) 0.017 (0.012) 57.968 (0.004)

1.63 (0.227) −0.346 (0.294) 1.073 (0.039)

1875.31 (320.63)

−1850.4 (345.37)

74.001 (0.075)

Plavink and Hurwitz, 1983 Gompertz Logistic Lopez Richards Sinusoidal

57.1 (9.71) 0.035 (0.015)

−51.8 (0.002)

2.17 (0.335) −0.125 (0.366) 1.21 0.041)

1820.9 187.1)

−1745.2 (55.9)

168.1 (0.059)

Hancock et al., 1995 (Ross 788) Gompertz Logistic Lopez Richards Sinusoidal

71.3 (3.51) 0.027 (0.003)

−49.4 (0.0004)

1.90 (0.056) −0.215 (0.072) 1.02 (0.063)

1852.8 (72.9)

1846.4 (11.8)

95.6 (0.017)

Govaerts et al., 2000 (Fast-H) Gompertz Logistic Lopez Richards Sinusoidal 1 2

64.8 (13.3) 0.034 (0.011) 40.72 (0.002)

1.92 (0.171) −0.135 (0.210) 0.9845 (0.022)

1228.73 (36.618)

-1213.6 (54.461)

74.01 (0.063)

W0 is initial weight, Wf is final weight, and a1 , a2, a3, b, K, and n are parameters. RMSE = Root mean square error, and R2 = Adjusted r squared.

and female strains of broiler chickens. Statistical analyses were performed using the non-linear procedure of MATLAB 7.13.0 employing the Nelder-Mead simplex algorithm.

RESULTS Limitations on space prevent the presentation of the results for all 28 time-course profiles. Therefore, as an example, model behavior for the 5 equations is illustrated in Figures 2 and 3 and the estimated parameters are given in Tables 3 and 4 for male and female data profiles, respectively. The results demonstrate that the

growth functions considered can be fitted to the data using the Nelder-Mead algorithm without difficulty or sensitivity to the choice of initial values for the parameters (Figures 2 and 3). The comparison of results (Tables 3 and 4) shows the sinusoidal function to be the most precise model in describing the growth data profiles from broilers used in this study. The values for the estimated parameters achieved with the logistic deviate the most from the values obtained with the other models. The logistic showed a trend to overestimate initial weights and underestimate final weights in both sexes for all data sources. The initial BW (W0 ) of the Gompertz, Richards, and

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

Lopez

82.2 (18.3) 144.7 (25.5) 62.2 (45.7) 50.3 (36.5) 55.01

233

SINUSOIDAL FUNCTION FOR GROWTH IN BROILERS Table 5. Comparison between the general goodness-of-fit of the models for given male (M) and female (F) data profiles based on various statistical criteria. Statistical criterion1

MLE MSE Errorparam

BIC

MLE MSE Errorparam AIC BIC

MLE MSE Errorparam AIC BIC

MLE MSE Errorparam AIC BIC

Gompertz

M F M F M F M F M F

− 126.1 − 114.2 5595.2 1884.4 74.8 43.4 266.3 242.3 257.1 233.2

M F M F M F M F M F

− 109.3 − 102.6 3254.1 1674.0 57.0 40.9 232.5 219.2 223.4 212.9

M F M F M F M F M F

− 348.7 − 339.6 526.3 407.0 22.9 20.2 709.4 691.2 708.2 690.0

M F M F M F M F M F

− 49.4 − 50.4 68.1 78.4 8.3 8.9 112.8 114.8 103.7 105.7

Logistic

Lopez

Richards

Sinus

Leeson and Summers, 1980 − 124.0 − 127.5 − 123.1 − 113.0 4613.5 6326.7 4243.1 1687.2 67.9 79.5 65.1 41.1 262.0 273.0 260.2 243.9 252.9 261.2 251.1 232.1

− 123.7 − 112.5 4487.1 1624.7 67.0 40.3 265.4 243.1 253.7 231.1

− 122.7 − 112.8 4086.5 1656.2 63.9 40.7 267.4 247.5 253.0 233.1

Plavink and Hurwitz, 1983 − 115.4 − 108.2 − 109.4 − 102.2 5997.4 2932.0 3306.0 1603.0 77.4 54.1 57.5 40.0 244.7 234.4 232.8 222.4 235.6 222.7 226.5 214.2

− 108.2 − 102.4 2939.9 1642.6 54.2 40.5 234.5 222.8 222.7 214.7

− 108.0 − 103.1 2870.9 1759.4 53.6 41.9 238.0 228.2 223.6 213.3

− 348.4 − 319.1 519.3 204.6 22.8 14.3 712.9 654.3 711.2 652.6

− 348.1 − 322.0 583.3 229.8 24.2 15.2 710.1 660.0 708.7 658.3

− 48.8 − 49.6 62.4 69.5 7.9 8.3 115.6 117.1 103.8 105.3

− 44.9 − 49.1 35.8 65.4 6.0 8.1 111.8 120.2 97.4 105.9

Hancock et al., 1995 − 422.8 − 363.4 − 418.6 − 322.0 5429.0 849.7 5196.4 225.6 73.7 29.1 72.1 15.0 857.6 740.9 849.2 660.0 856.3 739.4 847.9 658.3 Govaerts et al., 2000 (Fast-H) − 64.9 − 49.7 − 63.5 − 50.0 626.6 70.7 510.4 73.6 25.0 8.4 22.6 8.6 143.9 117.3 141.0 117.9 134.7 105.6 131.9 106.1

1 MLE = maximized log-likelihood, MSE = mean squared error, Errorparam = estimated parameters of the error variance, AIC = Akaike information criterion and BIC = Bayesian information criterion.

sinusoidal were lower than the logistic and Lopez estimates. The Gompertz, Richards, and sinusoidal W0 values were close to the expected initial average BW, while the logistic and Lopez initial weights were higher. For final (asymptotic) BW (Wf ), there were magnitude differences between the different models. A comparison between the functions showed lower estimated values of Wf for the logistic in most cases compared with the sinusoidal, Gompertz, Richards, and Lopez. Based on the results presented in Tables 3 and 4, female broilers reached their maximum BW gains earlier than the males in the same data source. A comparison among the models, based on the trend of the curves (Figures 2 and 3) showed that, except for the logistic, the growth functions examined gave a suitable fit to the profiles.

In general, a comparison between the models based on the estimated statistical criteria (Tables 5 and 6) indicated some relevant differences between the functions. The logistic equation gave poorer values of these criteria than the other growth functions. The chosen statistical criteria clearly demonstrate the suitability and superiority of the sinusoidal and Richards equations over the others. The sinusoidal equation was the best fitting model and the Richards the second best.

DISCUSSION Due to the nebulous nature of model development, modelers seem more willing to re-develop and adapt existing models for their own requirements than develop entirely new models. In summary, in this study a

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

AIC

Sex

234

DARMANI KUHI ET AL. Table 6. Comparison between the general goodness-of-fit of the models for the male (M) and female (F) data profiles based on various statistical criteria1 . Statistical 2

Criterion

MSE Errorparam AIC BIC

Sinus.3 vs. Gompertz

Sinus. vs. logistic

Sinus. vs. Lopez

Sinus. vs. Richards

M

F

M

F

M

F

M

F

86.7 73.3 46.7 86.7

86.7 80.0 53.3 73.3

100.0 100.0 86.7 93.3

86.7 93.3 100.0 100.0

86.7 60.0 86.7 86.7

66.7 80.0 46.7 53.3

86.7 53.3 53.3 86.7

53.3 53.3 40.0 40.0

Figure 4. Comparison between the general goodness-of-fit of the models to male data profiles using maximized log-likelihood estimates (MLE). A less negative value of MLE is an indication of better fit.

sinusoidal growth function was constructed and evaluated with regard to its ability to describe the relationship between live weight and age in broilers and compared with 4 standard growth functions, e.g., the Gompertz, logistic, Lopez, and Richards. The sinusoidal has a variable point of inflection, has the ability to approach final weight along either a gradual or an abrupt trajectory, and also has the ability to describe a wide range of hyperbolic shapes when there is no point of inflection. Comparison of the models based on their behavior (Figures 2 and 3) showed that, with exception of the logistic, the other growth functions gave a suitable fit to the data profiles. Here, the interesting choice lies between the sinusoidal and the Richards and Gompertz. Based on maximized log-likelihood, AIC and BIC criteria and depending on the sex (males better than females), the sinusoidal equation showed superiority to the other growth functions (Tables 5 and 6, and Figure 4). In many practical problems such as parameter estimation, the function values are uncertain or subject to noise. Therefore, a highly accurate solution is not necessary and in some instances may be impossible to compute. In such situations, all that is desired is an improvement in function value rather than a full optimization, which can converge to non-stationary points. To overcome these problems a global optimization method, the Nelder-Mead simplex, was used in this study to estimate the model parameters. This method does not

require any derivative information, which makes it suitable for problems with non-smooth functions. It has been widely used to solve parameter estimation and similar statistical problems where the function values are uncertain or subject to noise. It can also be used for problems with discontinuous functions, which occur frequently in statistics and experimental mathematics. The results of this study confirmed our assumption about the suitability of this optimization algorithm for curve fitting and growth parameter estimation.

REFERENCES Darmani Kuhi, H., E. Kebreab, and J. France. 2012. Application of the law of diminishing returns to partitioning metabolizable energy and crude protein intake between maintenance and growth in egg-type pullets. J. Appl. Poult. Res. 21:540–547. Darmani Kuhi, H., E. Kebreab, S. Lopez, and J. France. 2003. A comparative evaluation of functions for the analysis of growth in male broilers. J. Agric. Sci. 140:451–459. Darmani Kuhi, H., E. Kebreab, S. Lopez, and J. France. 2004. A comparative evaluation of functions for describing the relationship between live-weight gain and metabolizable energy intake in turkeys. J. Agric. Sci. 142:691–695. Darmani Kuhi, H., E. Kebreab, S. Lopez, and J. France. 2009. Application of the law of diminishing returns to estimate maintenance requirement for amino acids and their efficiency of utilization for accretion in young chicks. J. Agric. Sci. 147: 383–390. Darmani Kuhi, H., T. Porter, S. Lopez, E. Kebreab, A. B. Strathe, A. Dumas, J. Dijkstra, and J. France. 2010. A review of

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

1 Numbers in the table are the percentage of cases in which the fits of the sinusoidal model to data were superior to the other one according to the specified criteria in the rows. 2 MSE = mean squared error, Errorparam = estimated parameters of the error variance, AIC = Akaike information criterion and BIC = Bayesian information criterion. 3 Sinusoidal.

SINUSOIDAL FUNCTION FOR GROWTH IN BROILERS

Hooke, R., and T. A. Jeeves. 1961. Direct search solution of numerical and statistical problems. J. Assoc. Comput. Machin. 8:212– 229. Juretig, F., N. Buenos, and A. Argentina. 2012. Computing initial values for pharmacokinetic ML nonlinear regression via genetic algorithms and parallel genetic algorithms. Proceedings of the 2012 Pharmaceutical SAS Users Group Conference, Paper SP03. PharmaSUG, San Francisco, CA. http://www.pharmasug. org/proceedings/2012/SP/PharmaSUG-2012-SP03.pdf. Leeson, S., and J. D. Summers. 1980. Production and carcass characteristics of the broiler chicken. Poult. Sci. 59:786–798. Lister, D., T. Cowen, and R. A. McCance. 1966. Severe undernutrition in growing and adult animals. Br. Poult. Sci. 20:633– 639. Luenberger, D. G. 1984. Linear and Nonlinear Programming. 2nd Edition. Addison-Wesley, Reading, MA. MATLAB and Statistics Toolbox Release. 2012. The MathWorks, Inc., Natick, Massachusetts, United States. McCance, R. A. 1960. Severe under-nutrition in growing and adult animals. 1. Production and general effects. Br. Poult. Sci. 14:59– 73. Nelder, J. A., and R. Mead. 1965. A simplex method for function minimization. Comput. J. 7:308–313. Nocedal, J., and S. J. Wright. 1999. Numerical Optimization. Springer-Verlag, New York. Plavink, I., and S. Hurwitz. 1983. Organ weight and body composition in chicken as related to the energy and amino acid requirements: effects of strain, sex and age. Poult. Sci. 62:152–163. Sorensen, D. C. 1982. Newton’s method with a model trust region modification. SIAM J. Numer. Anal. 19:409–426. Thornley, J. H. M., and J. France. Mathematical Models in Agriculture, revised 2nd Edition. CAB International, Wallingford, UK; 2007; xvii+906pp. Torczon, V. 1993. On the convergence of pattern search algorithms. Tech. Rep. TR93-10, Rice University, Department of Mathematical Sciences, Houston, TX. Young, P. 1976. Optimization in the presence of noise. Pages 517– 573. In: Optimization in Action. L. C. W. Dixon (Ed.). Academic Press, London.

Downloaded from https://academic.oup.com/ps/article-abstract/97/1/227/4601754 by EKU Libraries user on 01 December 2019

mathematical functions for the analysis of growth in poultry. World. Poult. Sci. J. 66:227–239. Darmani Kuhi, H., F. Rezaee, A. Faridi, J. France, M. Mottaghitalab, and E. Kebreab. 2011. Application of the law of diminishing returns for partitioning metabolizable energy and crude protein intake between maintenance and growth in growing male and female broiler breeder pullets. J. Agric. Sci. 149:385–394. Dennis, J. E., and R. B. Schnabel. 1983. Numerical Methods for Unconstrained Optimization and Nonlinear Equations. PrenticeHall, Englewood Cliffs, NJ. Dumas, A, J. Dijkstra, and J. France. 2008. Mathematical modelling in animal nutrition: a centenary review. J. Agric. Sci. 146:123– 142. Faridi, A., S. L´ opez, H. Ammar, K. S. Salwa, A. Golian, J. H. M. Thornley, and J. France. 2015. Some novel growth functions and their application with reference to growth in ostrich. J. Anim. Sci. 93:2641–2652. Faridi, A., D. Murawska, A. Golian, M. Mottaghitalab, A. Gitoee, S. Lopez, and J. France. 2014. Alternative growth functions for predicting body, carcass, and breast weight in ducks: Lomolino equation and extreme value function. Poult. Sci. 93:1031–1042. Fitzhugh, H. A. 1976. Analysis of growth curves and strategies for altering their shape. J. Anim. Sci. 42:1036–1051. Fletcher, R. 1981. Practical Methods of Optimization. John Wiley & Sons, New York. France, J., M. S. Dhanoa, S. B. Cammell, M. Gill, D. E. Beever, and J. H. M. Thornley. 1989. On the use of response functions in energy balance analysis. J. Theor. Biol. 140:83–99. France, J., J. Dijkstra, and M. S. Dhanoa. 1996. Growth functions and their application in animal science. Ann. Zootech. 45:165– 174. Govaerts, T., G. Room, J. Buyse, M. Lippens, G. De Groote, and E. Decuypere. 2000. Early and temporary quantitative food restriction of broiler chickens. 2. Effects on allometric growth and growth hormone secretion. Br. Poult. Sci. 41:355–362. Hancock, C. E., G. D. Bradford, G. C. Emmans, and R. M. Gous. 1995. The evaluation of the growth parameters of six strains of commercial broiler chickens. Br. Poult. Sci. 36: 247–264.

235