El mejoramiento genético para incrementar la producción de leche en las regiones tropicales de climas cálidos se ha fundamentado, entre otras estrategias, en las cruzas de razas adaptadas (BI; Bos indicus) con razas especializadas de climas templados y fríos (BT; Bos taurus), aprovechando las ventajas del cruzamiento en diversos grados de encaste¹. En Cuba, a partir de los años 60 del siglo XX, se iniciaron cruzamientos entre Holstein y Cebú con el objetivo de obtener genotipos aptos para producir leche y adaptados a las condiciones tropicales de ambientes adversos; como resultado se generaron dos razas sintéticas^2,3, Siboney (5/8 Holstein y 3/8 Cebú) y Mambí (3/4 Holstein y 1/4 Cebú). Con respecto a Siboney, surgió como parte del Programa Nacional de Mejoramiento Genético para obtener genotipos menos exigentes que el Holstein, que desarrollaran un adecuado potencial productivo en condiciones tropicales con alimentación basada fundamentalmente en el pastoreo^4,5.

Para caracterizar la producción de leche en bovinos, la curva de lactancia (CL) se analiza a través de tres periodos, con cuatro componentes principales: 1) producción inicial; 2) fase ascendente o de incremento de la producción; 3) punto máximo o pico de producción; y, 4) tasa descendente o reducción de la producción, denominada persistencia^6,7,8. En la caracterización y análisis de la CL se han utilizado modelos no lineales (MNL), donde los coeficientes de regresión tienen una interpretación asociada a los estadios de la CL, o permiten derivar otros indicadores que describen la CL^9,10. El ajuste de los MNL se ha realizado con el análisis de regresión no lineal, utilizando métodos iterativos como Gauss, Marquardt y Newton, entre otros¹¹. Desde el punto de vista estadístico, la producción de leche con mediciones a través de la lactancia, representa un caso particular de medidas repetidas en el tiempo, generando una estructura de varianzas y covarianzas que puede ser analizada con modelos mixtos^12,13 con efectos positivos en los análisis para el ajuste de MNL. Los análisis de regresión no lineal, integrando efectos aleatorios relativos a los coeficientes de regresión, se realizan a través de método iterativo diferente (Firo, Gauss, Hardy e Isamp)^14,15 a los utilizados en los MNL con sólo efectos fijos.

El conocimiento de la CL permite orientar programas de manejo, alimentación y mejoramiento genético¹⁶. En sistemas de producción intensivos o estabulados, con razas especializadas como Holstein, el análisis de la CL se realiza a 305 días en leche; sin embargo, en los sistemas de producción tropicales, dado los efectos ambientales y las condiciones de manejo con repercusiones en la vida productiva de la vaca, la definición, análisis y caracterización de la lactancia se realiza a intervalos menores, como es el caso de 210 o 240 días^17,18,19. Para caracterizar y analizar la CL de bovinos Siboney, los objetivos del presente estudio fueron: a) seleccionar un MNL que se ajuste a la CL de bovinos Siboney; b) comparar los procedimientos de regresión no lineal, con MNL de solo efectos fijos vs MNL mixtos; c) estimar componentes o indicadores de la CL, con el modelo no lineal de mejor ajuste.

La base de datos fue proporcionada por el Centro de Investigaciones para el Mejoramiento Animal de la Ganadería Tropical, dependiente del Ministerio de Agricultura de la República de Cuba. Se analizaron 15,324 observaciones de producción diaria de leche (kg), en el intervalo de 1 a 240 días, de primeras lactancias en 2,809 vacas Siboney, con partos ocurridos de enero de 2000 a diciembre de 2012. Los datos provinieron de 28 hatos manejados en un sistema de pastoreo, ubicados entre los 20 y 23° N, y 74 y 85° O, con dos estaciones claramente definidas; la de lluvias de mayo a octubre, en la que ocurre del 70 al 80 % de la precipitación, y la estación seca de noviembre a abril; la temperatura media anual es de 23.1 °C y la humedad relativa de 60 a 70 % durante el día y de 80 a 90 % durante la noche²⁰. Para analizar y caracterizar la CL se evaluaron cinco MNL^7,9,21:

Donde: y_t= corresponde a la producción de leche (kg) en el día t (1 a 240); β₁, β₂ y β₃ = son los coeficientes de regresión que componen cada modelo; y ε_ij son los efectos aleatorios, asociados a los errores o residuales, de la j-ésima producción en el i-ésimo día de lactancia. Posteriormente, a partir de los β de cada modelo, se estimó la producción al inicio de la lactancia (PI), días para alcanzar la producción máxima (DPMAX), producción al pico o máxima (PMAX), y la producción total (PTOTAL) o acumulada a 240 días de lactancia, como indicadores de la CL^9,10,22.

Se realizaron dos análisis con el procedimiento NLMIXED de SAS, utilizando el método iterativo de Gauss^13,23: 1) MNL sin incluir efectos aleatorios relacionados con los coeficientes de regresión, sólo el efecto aleatorio (ε) de los residuales (ANA1); y, 2) MNL mixtos, se agregó el efecto aleatorio correspondiente al β₁ de cada modelo, más los residuales (ANA2). Los efectos aleatorios, no correlacionados, de ε_ij y β₁ se ajustan a las suposiciones de una distribución normal, con media igual a cero y varianzas de σ²_e y σ²β₁, respectivamente; no se incluyeron los efectos aleatorios de los coeficientes de regresión β₂ y β₃, dado que sus componentes de varianza son de muy baja magnitud, y tienden a cero. La selección del modelo con mejor ajuste se realizó en función de seis criterios^24,25,26:

Donde: pmi= producción medida en el i-ésimo día de lactancia; pei= producción estimada para el i-ésimo día de lactancia; n= número total de observaciones; sce= suma de cuadrados de los residuales; sct= suma de cuadrados total; Log lik= logaritmo de la función de verosimilitud; k= número de parámetros en el modelo. El EPP analiza la relación que existe entre la producción medida y la producción estimada para el i-ésimo día de lactancia, y en función del signo, el MNL sobreestima (+) o subestima (-) las predicciones. Para EPP, VEP, AIC y BIC, el modelo con el menor valor se consideró como el de mejor ajuste; a diferencia del R², el modelo con el valor más alto representa el mejor ajuste. El estadístico DW analiza las auto correlaciones en los errores, con tres planteamientos: a) si 2<DW<4 existe auto correlación negativa; b) si 0<DW<2 indica ausencia de auto correlación; y, c) si DW≤0 indica que existe auto correlación positiva.

En el Cuadro 1, se presentan los resultados de los criterios estadísticos para la selección del modelo de mejor ajuste; para EPP, VEP y DW los resultados son similares a través de los dos análisis realizados. Con base en el EPP y el estadístico DW, las predicciones de todos los MNL tendieron a subestimar la información analizada, con ausencia de auto correlación entre los residuales. Los MNL mixtos mostraron el mejor ajuste con base en tres resultados del Cuadro 1: a) las varianzas residuales de los MNL del ANA2, fueron menores (hasta 82 % para los modelos de SIK y WOD) a las estimadas en los ANA1; b) por consiguiente, los R² del ANA2 (96.7 % en promedio) fueron superiores a los obtenidos en el ANA1 (87 % en todos), en un intervalo de 9 a 12 %; c) y conjuntamente, los criterios AIC y BIC, mostraron mejor ajuste en los MNL del ANA2; posteriormente, dentro del ANA2, el modelo de mejor ajuste fue WOD, seguido de SIK y BRO, en función de los criterios AIC, BIC y R². Las estimaciones de PI, DPMAX, PMAX y PTOTAL se realizaron con los resultados del ANA2 (Cuadro 1). Diferentes autores coinciden en publicar sobre el MNL de WOD, como de mejor ajuste estadístico a través de diversos criterios de selección o jerarquización, o en su caso más popular, para describir y analizar las CL de bovinos para leche; los coeficientes de regresión que lo conforman se asocian a diferentes componentes de la CL, y permiten derivar otros indicadores de la producción, parcial o total, a través del tiempo^7,8,27; García-Muñiz et aX²¹⁾, reportaron que el MNL de WOD resultó el idóneo, a través de 16 MNL evaluados, en CL de seis genotipos de bovinos cruzados, con diferentes proporciones de BT y BI; así mismo, en bovinos para carne, el modelo de WOD ha resultado apto, con ajuste idóneo, para modelar la producción de leche a través de diferentes esquemas de muestreo^28,29.

Los β₁, β₂ y β₃ de los MNL están asociados a la producción al inicio de la lactancia, a la fase ascendente hasta la producción máxima, así como a la tasa descendente o persistencia, respectivamente. Dada la combinación de signos en β₂ y β₃, los modelos de WOD y WIL permiten ajustar cuatro tipos de curvas unimodales³⁰: 1) curva estándar o típica: WOD + β₂,- β₃ y WIL - β₂,- β₃; 2) curva continuamente decreciente (atípica): WOD - β₂,- β₃ y WIL + β₂,- β₃; 3) curva estándar invertida (U): WOD - β₂,+ β₃ y WIL + β₂,+ β₃; y, 4) continuamente incrementando: WOD + β₂,+ β₃ y WIL -β₂,+ β₃. En la Figura 1 se presenta la CL estándar o típica con base en el modelo de WIL, y la CL continuamente incrementando con el modelo de WOD. Las estimaciones de los componentes de la CL, a través de los MNL evaluados, fueron similares para PMAX y PTOTAL, y con marcadas diferencias para PI y DPMAX (Cuadro 1). Los modelos de COB y BRO, expresaron valores subestimados para puntos iníciales de la lactancia, como es el caso de PI y DPMAX (1.8 kg y 20 días, respectivamente), en contraste a las altas estimaciones de β₁; lo anterior puede señalar ciertos problemas de estos modelos para ajustar la parte inicial de la CL.

Figura 1 Descripción de la curva de lactancia de bovinos Siboney, ajustada a 240 días con base en los modelos de Wiltmink (a; WIL) y Wood (b; WOD) 

Los valores promedios para PMAX y PTOTAL fueron 7.8 y 1,660 kg, respectivamente; estos resultados son superiores a los previamente reportados^5,31,32 en bovinos Siboney; además de los reportados por Rodríguez y Ponce de León¹⁹ en otros genotipos producto de las cruzas de BI con BT; sin embargo, son inferiores a los reportados en bovinos Mambí^33,34.

Los sistemas de producción de bovinos lecheros en condiciones tropicales cálidas dependen de factores relacionados con el animal, ambiente, alimentación y tecnología de producción. En el contexto del animal, la variabilidad en los niveles productivos se puede atribuir al genotipo, en función de la proporción de genes BI y BT que constituyen al individuo^35,36; en el contexto de tecnología, la producción de leche es afectada por la interacción del becerro con la vaca y por el manejo del amamantamiento de la cría¹⁸. La alimentación en los sistemas tropicales se fundamenta en pastoreo de gramíneas forrajeras, y en algunas ocasiones se complementa con leguminosas forrajeras³⁷; sin embargo, los efectos del ambiente se caracterizan por la prevalencia de parásitos internos y externos, temperaturas y niveles de humedad altos que repercuten en estrés del animal; así como en la calidad y disponibilidad del forraje^1,38. En bovinos Siboney, el modelo de Wood fue el que mejor se ajustó para caracterizar la curva de lactancia a 240 días; presentando una curva continuamente incrementando, con una producción máxima de 7.71 kg a los 40 días, así como una producción total acumulada de 1,653 kg. Los modelos no lineales mixtos presentaron los mejores ajustes, en función de los criterios de información Akaike y Bayesiano; además de reducir la varianza residual, con aumentos del 11 % en el coeficiente de determinación.