Proposta de análise genética de curvas de crescimento de bovinos por meio do algoritmo SAEM

Proposal for genetic analysis of growth curves of cattle through the SAEM algorithm

N.A.M. Silva^I; Â.M.Q. Lana^II; F. Fonseca e Silva^III; R.R. Lima^IV; M.A. Silva^II; J.A.G. Bergmann^II

^IUniversidade Federal de Uberlândia - Uberlândia, MG
^IIUniversidade Federal de Minas Gerais - Belo Horizonte, MG ]]> ^IIIUniversidade Federal de Viçosa - Viçosa, MG
^IVUniversidade Federal de Lavras - Lavras, MG

RESUMO

Compararam-se duas diferentes metodologias na avaliação genética de curvas de crescimento de animais Nelore: o algoritmo SAEM e o método Two-step. Para a implementação dessas metodologias, foram utilizados o modelo de crescimento de Brody modificado e o modelo touro. A diferença entre o SAEM e o Two-step é que o algoritmo SAEM estima simultaneamente parâmetros do modelo e efeitos genéticos e ambientais, e o método Two-step faz esse processo de estimação em duas etapas distintas. Mais ainda, o algoritmo SAEM utiliza o método de máxima verossimilhança, e o do Two-step o de máxima verossimilhança restrita. Foram obtidos, com base nas metodologias testadas, além das estimativas de efeitos fixos e parâmetros genéticos, os valores genéticos preditos para os touros avaliados. A partir dos valores genéticos preditos, foram obtidas as curvas genéticas para os touros. O algoritmo SAEM mostrou-se consistente na estimação dos efeitos fixos e na predição dos efeitos aleatórios, apresentando-se como uma alternativa viável para avaliação genética de animais Nelore.

Palavras-chave: componentes de (co)variância, algoritmo estocástico EM, curvas genéticas

ABSTRACT

Two methodologies in genetic evaluation of growth curves of Nellore cattle were compared: the SAEM algorithm and the Two Step method. To implement these methodologies the Brody modified growth curve and the sire model were used. The difference between the SAEM and the Two Step is that SAEM estimates jointly the parameters of the model and genetics and environmental effects and the Two Step method does this process in two independent steps. Estimates of the fixed effects and genetics parameters, and prediction breeding values for the sires were obtained from the methodologies. From the breeding values genetic curves were obtained for the sires. The SAEM algorithm proved consistent in the estimation of fixed effects and prediction of random effects.

Keywords: (co)variance components, stochastic algorithm EM, genetics curves

INTRODUÇÃO

Modelos de regressão não linear têm se mostrado adequados para descrever curvas de crescimento de animais, pois apresentam parâmetros que podem ser interpretados biologicamente (Silva et al., 2010). Entre os modelos utilizados para descrever o crescimento, o modelo de Brody (1945) apresenta-se adequado para descrever o crescimento de bovinos Nelore. Este modelo, segundo Forni (2007), é o mais utilizado para descrever o crescimento de bovinos de corte, pela facilidade de convergência dos parâmetros, pelo menor esforço computacional e pela qualidade do ajuste.

As estimativas dos parâmetros dos modelos de crescimento em programas de seleção, normalmente, são obtidas por meio de uma metodologia frequentista (método Two-step), que considera duas fases distintas (Garnero et al., 2005). De acordo com Blasco et al. (2003), a estimação em duas fases distintas não é apropriada, pois a influência de efeitos genéticos e ambientais não é considerada no processo de estimação dos parâmetros de crescimento, o que pode caracterizar em perda de informação. De acordo com estes autores, uma forma eficiente de se realizar esta análise é por meio da inferência bayesiana.

Jaffrézic et al. (2006) relataram que a inferência bayesiana realmente apresenta grande eficiência na estimação simultânea de parâmetros de crescimento e de parâmetros genéticos, porém é muito complexa, uma vez que exige a obtenção de distribuições condicionais completas e aplicação de algoritmos, cujas convergências são de difícil constatação. Como alternativa de estimação simultânea, estes autores sugerem uma extensão do algoritmo EM com o uso de um processo estocástico, denominado algoritmo Stochastic Approximation EM algorithm (SAEM), para se obter estimativas pelo método da máxima verossimilhança, ou seja, recomendam utilizar o enfoque frequentista. Ainda, segundo esses autores, o SAEM é bastante robusto para valores iniciais e tem muita rapidez na convergência das estimativas de máxima verossimilhança.

Em melhoramento genético de bovinos de corte, conhecer somente o comportamento genético dos parâmetros da curva não é uma alternativa que esteja sendo implementada na prática em programas de melhoramento (Albuquerque e Fries, 1998). Assim, de acordo com estes autores, além de conhecimento do perfil genético de tais parâmetros, é importante conhecer a denominada curva de crescimento genética dos animais, a qual permite predizer o valor genético dos animais em qualquer ponto do intervalo de tempo em que as medidas de peso e a idade foram tomadas, possibilitando, dessa forma, realizar a seleção em qualquer idade desejada.

Os objetivos da realização deste trabalho foram: obter, a partir da utilização do algoritmo SAEM e do método Two-step, os componentes de variância dos parâmetros com interpretação biológica do modelo não linear de crescimento de Brody modificado, bem como estimar os valores de herdabilidade de cada parâmetro do modelo e a correlação genética entre os parâmetros. Além disso, serão obtidas as curvas genéticas individuais, relacionadas a cada uma das metodologias empregadas, a partir dos valores genéticos preditos para os parâmetros do modelo de crescimento utilizado.

MATERIAL E MÉTODOS

O modelo de crescimento não linear utilizado para avaliar o crescimento dos animais, foi o de Brody. Esse modelo é descrito da seguinte forma (Jaffrézic et al., 2006): Y_ij = A_i - B_ie - K_it_ij + ε_ij, em que: y_ij é o peso do indivíduo i tomado no tempo j, A_i, B_i, K_i são os parâmetros específicos de cada indivíduo i, sendo que Ai representa a estimativa de peso adulto, K_i a taxa de maturidade e B_i não possui interpretação biológica. Como sugerido por Jaffrézic et al. (2006), pelo fato de B_i não possuir interpretação biológica, uma reparametrização foi utilizada, de tal forma que B_i = Ai - y_i1, em que y_i1 representa o peso observado do animal i na idade 1, ou seja, o peso ao nascimento do animal. Outros autores, Tedeschi et al. (2000) e Silva et al. (2002), já usaram essa reparametrização, chamando esse modelo de Brody modificado: como tij é a idade correspondente à observação, sendo t_ij= t_ij/100000, é o erro aleatório, independentemente e normalmente distribuído, com ε_ij~N(0, σ²e).

O modelo touro prediz o efeito do pai sobre o desempenho médio de sua progênie, desconsiderando o efeito da mãe. Neste modelo, os parâmetros A e K de cada indivíduo, filhos de cada touro, são determinados por um efeito fixo β (média do parâmetro, já que não foram considerados grupos contemporâneos) e por um efeito genético u atribuído ao touro, além de um efeito aleatório não controlado denotado por e. Considerando-se essas definições, o vetor Φ_i para o individuo i é decomposto como: Φ_i = X_iβ + Z_iu + e_i (i= 1, ..., N), em que: Φ_i é um vetor de parâmetros 2Nx1; Φ_i =[A1 ,A2 , ... , AN ,K1, K2,... , KN ]'; X _i é uma matriz de incidência dos efeitos fixos (2Nx2); β é um vetor de efeitos fixos (2x1); Z_i é uma matriz de incidência dos touros (2NxS), sendo S o número de touros; u é um vetor de efeitos genéticos atribuídos aos touros (S x 1); u_i ~ N(0,⊗G), sendo A a matriz de parentesco entre os touros e G (2x2) a matriz de (co)variância genética entre os parâmetros da curva. No modelo touro, os us preditos equivalem à metade do valor genético do indivíduo (touro); e_i é um vetor de resíduos referentes ao modelo (2Nx1); e_i~N(0,⊗I), sendo P (2x2) uma matriz de covariância residual entre os parâmetros do indivíduo i e I uma matriz de identidade NxN.

Para a obtenção dos componentes de variância provenientes do ajuste do modelo touro, os dados de bovinos Nelore foram analisados por meio de duas metodologias: SAEM e Two-step.

De acordo com Jaffrézic et al. (2006), a ideia geral do algoritmo SAEM é substituir a fase de esperança do algoritmo EM por uma aproximação estocástica. Esta substituição deve ser feita porque o cálculo da esperança condicional da função de log-verossimilhança não possui uma forma analítica de ser calculada (Panhard e Samson, 2008). Ainda segundo Jaffrézic et al. (2006), na estrutura do algoritmo EM, para a utilização do SAEM, na análise genética de curvas de crescimento, o vetor de parâmetro θ a ser estimado será θ =(β,G,P,σ²_e). Os dados faltantes serão Φ e u, e o vetor de dados completos é então descrito por x=(y, z)', sendo y o vetor de dados observados e z = (Φ, u).

Cada iteração do algoritmo SAEM poderá ser sumarizada com os seguintes passos:

1. simulação - foram simuladas 10 cadeias aleatórias (L=10) de Φ e u. Para simulação de cada cadeia foi utilizado o algoritmo Metropolis Hasting (M). As cadeias simuladas são substituídas iterativamente nas aproximações estocásticas dadas no passo 2;

2. aproximação estocástica - é a substituição da esperança condicional Q(θ|θ^k) = E[logL(θ;y,z)|y,θ^k] pela seguinte aproximação estocástica:

Segundo Jaffrézic et al. (2006), para cada parâmetro a ser estimado, β, G, P, e σ²_e σ²_e, têm-se as seguintes aproximações estocásticas, respectivamente:

3. maximização de Q(θ|θ(k)) - os parâmetros que serão estimados, β, G, P, e σ2, serão maximizados com o uso das suas respectivas aproximações estocásticas dadas no passo 2. Assim, as estimativas desses parâmetros serão obtidas em cada iteração da seguinte forma, de acordo com Jaffrézic et al., 2006:

,em que: N é o número de indivíduos com observações; N_tot é o número total de observações; e N_a é o número de animais na matriz de parentescos.

Após a obtenção dos valores que maximizam a função de máxima verossimilhança neste último passo, estes são atualizados no passo 1, dando origem ao processo iterativo. Tal processo é então repetido até se obter a convergência.

O método frequentista Two-step é o mais utilizado para obtenção dos componentes de variância. Vários autores utilizaram esse método em estudos com bovinos: Silva et al. (2002), Talhari et al. (2003), Garnero et al. (2005), Santoro et al. (2005) e Silva (2009),. Divide-se em duas etapas independentes para obtenção dos componentes da variância. A primeira etapa consiste na obtenção dos parâmetros do modelo de crescimento. Na segunda etapa, as estimativas obtidas para os parâmetros do modelo da etapa 1 são utilizadas como variáveis dependentes em um modelo touro contendo os efeitos fixos e aleatórios. A partir desse modelo touro, podem-se obter os componentes de variâncias necessários aos cálculos das herdabilidades e das correlações.

A curva de crescimento genética dos touros foi obtida de acordo com Pong-Wong e Hadjipavlou (2009), que utilizaram o modelo de regressão não linear de Gompertz, para fazer a avaliação genômica em ovelhas, da seguinte forma: , em que: é o vetor de valores genéticos acrescidos da média para os pesos do touro i em cada idade j; Ua* e Uk* são os vetores de valores genéticos acrescidos da média para as estimativas dos parâmetros A e K (Ua* = µ_a + Ua e Uk* = µ_k + Uk ); Ua são valores de us relativos ao parâmetro A; Uk são valores de us relativos ao parâmetro K; PN é a média do peso ao nascimento para os filhos touro i; tij são as pesagens simuladas j relacionadas ao touro i, sendo que t variou de 0 a 550 dias de idade, com um intervalo de 10 em 10 dias.

Além da estimação da curva genética, foi estimada também uma curva média para os touros, que pode ser descrita, então, da seguinte forma: , em que µ_a e µ_k são os efeitos fixos para A e K, respectivamente.

A partir da diferença entre os valores gerados pela curva média ((y)) e a curva genética estimada para cada touro (y), têm-se os valores genéticos preditos para cada touro em qualquer idade de interesse (do nascimento aos 550 dias). Por meio dessa técnica, foram obtidos os valores genéticos desses touros em idades-padrão: nascimento (0 dias), desmama (210 dias), ano (360 dias) e sobreano (540 dias). Os valores genéticos obtidos para essas idades foram correlacionados entre as diferentes metodologias testadas a fim de se verificar a semelhança delas na predição de tais valores.

Além dos valores genéticos obtidos a partir da metodologia descrita, também foram obtidos valores genéticos por meio do ajuste do modelo touro em análises unicaracterísticas para os pesos observados nas idades próximas às idades-padrão. Tais valores genéticos foram correlacionados pelo coeficiente de Spearman com os valores genéticos obtidos pelas curvas genéticas para as duas metodologias testadas.

A forma de aplicação dos métodos SAEM e Two-step se deu de quatro maneiras distintas, conforme descrição a seguir: 1a análise) SAEM1 - os dados dos 410 animais Nelore foram analisados a partir do algoritmo SAEM. Pelo fato de esse algoritmo realizar uma análise conjunta, não há problemas de convergência dos parâmetros do modelo, e todos os animais foram incluídos no processo de estimação; 2a análise) Two- step 1 (TS1) - dos 410 animais utilizados na primeira etapa de ajuste deste método, apenas 326 obtiveram convergência para os parâmetros do modelo. Assim, para a segunda etapa do processo, somente foram utilizadas as estimativas de parâmetros de 326 animais, reduzindo, dessa forma, o número de pais para avaliação no modelo touro para 89; 3a análise) SAEM2 - mediante a diminuição no número amostral no TS1, e para fins de comparação das metodologias utilizando o mesmo número de animais, o algoritmo SAEM foi empregado para analisar apenas os dados dos mesmos 326 animais que convergiram no TS1; 4a análise) Two-step 2 (TS2) - ainda sob o enfoque de comparar as metodologias utilizando os mesmos números de animais, e com o intuito também de verificar se a forma de estimação conjunta dos parâmetros do modelo influencia na estimação dos componentes de variância, uma segunda alternativa para o Two-step foi avaliada. Nesta análise, os parâmetros do modelo (A e K) estimados pelo SAEM 1 foram utilizados como as variáveis dependentes no modelo touro, de modo que, a partir destes parâmetros, já estimados conjuntamente pelo SAEM 1, fosse possível obter os componentes de variância.

RESULTADOS E DISCUSSÃO

A partir dos resultados apresentados na Tab. 1, percebe-se que, em relação aos efeitos fixos estimados, não houve grande diferença entre as análises realizadas. Os valores encontrados ficaram condizentes com Silva et al. (2002), que também utilizaram a função de Brody modificada com bovinos Nelore. Quanto aos componentes de variância, pode-se perceber que o SAEM1 apresentou diferença em relação ao TS1, principalmente quanto aos componentes residuais. Isso se deve provavelmente ao fato de TS1 ser executado em duas etapas distintas, assim os erros provenientes do ajuste do modelo de crescimento na primeira etapa não são considerados na estimação dos componentes de variância na segunda etapa, o que pode também estar elevando essa variação residual. Como no SAEM1 a estimação é simultânea, há uma correção nas estimativas dos componentes de variância, diminuindo a variação residual. O mesmo comportamento é observado quando se compara TS1 e SAEM2. Quando se compara SAEM1 com TS2, as estimativas dos componentes de variância obtidas se aproximam. Esse fato elucida que a diferença na estimação dos componentes de variância é realmente proveniente do erro de ajuste da primeira etapa, pois o TS2 utiliza as estimativas de A e K do SAEM1 para se estimar os componentes de variância.

É importante também considerar que outro fator que pode estar elevando as diferenças entre os componentes de variância estimados é o fato de o SAEM e o TS considerarem metodologias diferentes para estimação desses componentes. O SAEM utiliza o método de ML, e o TS o de REML. Por esse motivo, as estimativas do SAEM podem estar sendo subestimadas, pois este já é um comportamento esperado quando se utiliza a ML. Meza et al. (2007) comprovam esse fato em um estudo de simulação de dados, em que estes autores compararam estimativas de componentes de (co)variância dos parâmetros da curva de crescimento de Gompertz pelos métodos ML e REML. Esses autores obtiveram, conforme o esperado, menores valores quando utilizada a ML, principalmente para os componentes genéticos em relação aos residuais. Assim, pode-se concluir que as diferenças observadas podem, na verdade, ser menores que as observadas na Tab. 1, pelo viés da estimativa, quando se considera o uso da ML.

Ainda, pela Tab. 1, tendo em vista as estimativas de herdabilidade obtidas para o parâmetro A, pode-se dizer que os valores são concordantes com as discussões relativas às semelhanças entre as análises SAEM1 e TS2. Silva et al. (2002), em revisão sobre valores de herdabilidade de parâmetros de modelos de crescimento aplicados a bovinos de corte, relataram valores de herdabilidade de A variando de 0,04 a 0,46, e de K variando de 0,08 a 0,41. Os valores de correlação genética, com exceção do método TS2, também condizem com o encontrado por esses autores, sendo altos e negativos. Em relação ao valor encontrado para esse parâmetro na análise TS2, de acordo com Silva et al. (2010), é um comportamento de difícil explicação e, possivelmente, deve-se a problemas de amostragem.

De forma geral, pode-se perceber, pela Tab. 2, que, embora tenham sido observadas diferenças nos componentes de variância e nos parâmetros genéticos (Tab. 1), os valores genéticos para o parâmetro A tenderam a ser semelhantes no ranking de seus valores, uma vez que as correlações de Spearman variaram de 0,63 a 0,82, e foram todas significativas. Percebe-se, ainda, que as análises que tiveram maiores semelhanças na estimação dos componentes de variância (Tab. 1) foram as que apresentaram maiores valores de correlação entre os valores genéticos para A, confirmando, então, a concordância do SAEM em seu processo de estimação e a baixa precisão do TS na primeira etapa do método. Para os valores genéticos do parâmetro K, pode-se observar que eles também foram significativos, embora de menor magnitude, variando de 0,51 a 0,87.

Pelas Fig. 1 e 2, nota-se, de forma geral, que as curvas genéticas obtidas pelas análises TS1 e TS2 apresentam diferenças em relação à variação das curvas em torno da curva média, e que essa diferença marcante não é observada nos resultados obtidos pelas análises SAEM1 e SAEM2, corroborando as discussões nas quais se verifica a maior estabilidade do método SAEM.

Foram obtidos para cada touro os valores genéticos em idades-padrão, sendo que tais valores representaram as diferenças entre os valores genéticos da curva genética individual e os valores médios. Esses foram correlacionados por meio do coeficiente de Spearman, para todas as análises realizadas, como mostra a Tab. 3. Para todas as idades avaliadas, as maiores correlações foram encontradas entre SAEM1 e SAEM2 e por SAEM1 e TS2. Tais resultados também confirmam, mais uma vez, a semelhança, ou a já referida estabilidade, das análises baseadas no algoritmo SAEM.

Os coeficientes de correlação (Tab. 3) tenderam a se tornar maiores à medida que as idades-padrão aumentaram. Dessa forma, pode-se considerar que a seleção feita em idades mais avançadas - a partir de um ano - tendeu a ser menos influenciada pelo método em questão. Além disso, pode-se relacionar tal resultado com as curvas mostradas nas Fig. 1 e 2, nas quais se notou que a variância dos valores genéticos em cada tempo aumentou com o avançar da idade.

CONCLUSÕES

O algoritmo SAEM mostrou-se adequado para a proposta de avaliação genética de bovinos, pois apresentou estimativas mais estáveis quando comparadas ao TS, além de ter predito valores genéticos para idades-padrão correlacionados significativamente com aqueles obtidos via análises unicaracterísticas usuais. É importante considerar que este trabalho serviu apenas para ilustrar a metodologia do algoritmo SAEM, aqui descrita como sendo um método alternativo de avaliação genética de curvas de crescimento.

REFERÊNCIAS

ALBUQUERQUE, L.G.; FRIES, L.A. Modelos de regressão aleatória: possível aplicação na seleção de animais precoces. In: CONGRESSO BRASILEIRO DAS RAÇAS ZEBUÍNAS, ABCZ e SEBRAE, 3., Uberaba. Anais... Uberaba, 1998. p.129-144         [ Links ]

BLASCO, A.; PILES M.; VARONA, L. A Bayesian analysis of the effect of selection for growth rate on growth curves in rabbits. Genet. Select. Evol., v.35, p.21-41, 2003.         [ Links ]

BRODY, S. Bioenergetics and Growth. New York: Rheinhold Publishing, 1945. 1023p.         [ Links ]

FORNI, S. Analise da curva de crescimento de bovinos da raça nelore utilizando funções não lineares em análises bayesianas. 2007. 75f. Tese (Doutorado em Zootecnia) - Universidade Estadual Paulista, Jaboticabal, SP.         [ Links ]

GARNERO, A.V.; MARCONDES, C.R.; BEZERRA, L.A.F. et al. Parâmetros genéticos da taxa de maturação e do peso assintótico de fêmeas da raça Nelore. Arq. Bras. Med. Vet. Zootec., v.57, p.652- 662, 2005.         [ Links ]

JAFFRÉZIC, F.; MEZA, C.; LAVIELLE, M.; FOULLEY, J.L. Genetic analysis of growth curve using the SAEM algorithm. Genet. Sel. Evol., v.38, p.583-600, 2006         [ Links ]

MEYER, K. WOMBAT - Digging deep for quantitative genetic analyses by restricted maximum likelihood . In: WORLD CONGRESS ON GENETICS APPLIED TO LIVESTOCK PRODUCTION, 8., 2006, Belo Horizonte. Proceedings... Belo Horizonte. 2006. CD-ROM.         [ Links ]

MEZA, C.; JAFFRÉZIC, F.; FOULLEY, J.L. REML Estimation of Variance Parameters in Nonlinear Mixed Effects Models Using the SAEM Algorithm. Biometrical J., v.49, p.876-888, 2007.         [ Links ]

PANHARD, X.; SAMSON, A. Extension of the SAEM algorithm for nonlinear mixed models with 2 levels of random effects. Biostatistics, v.10, p.121-135, 2009.         [ Links ]

PONG-WONG, R.; HADJIPAVLOU, G. A two-step approach combining the Gompertz growth with genomic selection for longitudinal data. BMC proceedings, v.4, suppl. 1:S4, p.4, 2010.         [ Links ]

STATISTICAL Analysis System user's guide. SAS Institute Inc. Version 9.1.ed. Cary: SAS Institute, USA, 2003.         [ Links ]

SANTORO, R.S.; BARBOSA, S.B.P.; BRASIL, L.H.A.; SANTOS, E.S. Estimativas de parâmetros de curvas de crescimento de bovinos zebu, criados no estado de Pernambuco. Rev. Bras. Zootec., v.34, p.2262-2279, 2005.         [ Links ]

SILVA, F.L. Curvas de crescimento e produtividade de vacas nelore e cruzadas, de diferentes tipos biológicos, em sistema de produção intensiva. 2009. 89p. Dissertação (Mestrado em Zootecnia) - Escola Superior de Agricultura " Luiz de Queiroz", Piracicaba, SP.         [ Links ]

SILVA, N.A.M.; LIMA, R.R.; SILVA, F.F.E. et al. Modelo Hierárquico Bayesiano aplicado na avaliação genética de curvas de crescimento de bovinos de corte. Arq. Bras. Med. Vet. Zootec., v.62, p.647-654, 2010         [ Links ]

SILVA, F.F.; AQUINO, L.H.; OLIVEIRA, A.I.G. Estimativas de parâmetros genéticos de curva de crescimento de gado Nelore (Bos indicus). Cienc. Agrotecnol., v.26, p.1562-1567, 2002.         [ Links ]

TALHARI, F.M.; ALENCAR, M.M.; MASCIOLI, A.S. et al. Correlações genéticas entre características produtivas de fêmeas em um rebanho da raça Canchim Rev. Bras. Zootec., v.32, p.880-886, 2003.         [ Links ]

TEDESCHI, L.O.; BOIN, C.; NARDON, R.F.; LEME, P.F. Estudo da curva de crescimento de animais da raça guzerá e seus cruzamentos alimentados a pasto, com e sem suplementação. 1. análise e seleção das funções não-lineares. Rev. Bras. Zootec., v.29, p.630-637, 2000.         [ Links ]

Recebido em 27 de outubro de 2011 ]]> Aceito em 13 de março de 2012

Email: natascha@famev.ufu.br