Analisis del articulo: Phylogeny and metabolic scaling in mammals
Maria Camila Chaves y Sebastian Hidalgo
2026-05-29
1 Abstract
The relationship between metabolism and body mass in mammals has been highly debated throughout allometry studies. Traditional studies propose that the universal allometric exponent corresponds to 2/3 or 3/4, however in the article Phylogeny and Metabolic Scaling in Mammals by Isabella Capellini et al. (2010) proposes something different; the calculation of this relationship must incorporate the phylogenetic history of the species to obtain more precise estimates. In this study, the methodology of the original article was recreated by analyzing basal metabolism (BMR) and field metabolism (FMR) data through linear regressions (log-log), as well as ordinary least squares (OLS) models, independent phylogenetic contrasts (PIC) and PGLSk models (which they propose is more optimal for this analysis. All the analyzes carried out corresponded to linear regression models; however, each method incorporated the phylogenetic dependence between each other in a different way. species. While the OLS models assumed total statistical independence, the PIC and PGLS models explicitly incorporated the evolutionary history of the species through phylogenetic contrasts and variance-covariance matrices derived from the phylogeny obtained from the Bininda-Emonds super-tree of mammals. The purpose of this document is to review the methodology and execution of the article, with the intention of seeing how close the results can be, finally confirming the importance of incorporating phylogeny in studies of species. metabolic scaling and supports the use of PGLSk models as a more robust approach compared to traditional methods.
2 Abstracto
La relación entre el metabolismo y masa corporal en los mamíferos ha sido altamente debatida a lo largo de los estudios de la alometria. Los estudios tradicionales proponen que el exponente alometrico universal corresponde a 2/3 o 3/4, no obstante en el artículo Phylogeny and Metabolic Scaling in Mammals de Isabella Capellini et al. (2010) se propone algo distinto;el cálculo de esta relación debe incorporar la historia filogenética de las especies para obtener estimaciones más precisas. En este estudio,se recreó la metodología del artículo original mediante el análisis de datos de metabolismo basal (BMR) y metabolismo en campo (FMR) por medio de regresiones lineales (log-log), así como modelos de mínimos cuadrados ordinarios (OLS), contrastes filogeneticos independientes (PIC) y modelos PGLSk (el cual plantean que es mas optimo para este análisis. Todos los análisis realizados correspondieron a modelos de regresión lineal; sin embargo, cada método incorporó de manera diferente la dependencia filogenética entre especies. Mientras los modelos OLS asumieron independencia estadística total, los modelos PIC y PGLS incorporaron explícitamente la historia evolutiva de las especies mediante contrastes filogenéticos y matrices de varianza-covarianza derivadas de la filogenia obtenida del super-tree de mamíferos de Bininda-Emonds. La finalidad de este documento es hacer una revisión de la metodología y ejecución del articulo, con la intención de ver que tan cercanos pueden ser los resultados, finalmente confirmar la importancia de incorporar la filogenia en estudios de escalamiento metabólico y respalda el uso de modelos PGLSk como una aproximación más robusta a comparación de los métodos tradicionales.
3 Introducción
Contexto biológico/ecológico y pregunta de investigación del artículo original. La alometria biológica, es entre muchas cosas esa relación entre el tamaño corporal del del animal con su tasa metabólica.Esta relación se sabe que no aumenta de manera proporcional al tamaño corporal, por lo que se usan ecuaciones de escalamiento para su descripción.Esta relación ha sido considerada fundamental para comprender múltiples procesos ecológicos y evolutivos, debido a que el metabolismo influye directamente sobre variables como el crecimiento, la reproducción, el uso de energía, la dinámica poblacional y las estrategias de vida de los organismos.
Antes de este articulo habían dos exponentes alometricos que estaban dominando el debate, el primero corresponde al modelo de superficie-volumen, el cual propone que el metabolismo escala con la masa corporal siguiendo un exponente cercano a 2/3, debido a las limitaciones impuestas por la disipación de calor a través de la superficie corporal. Por otro lado, la conocida “ley de Kleiber” propone un exponente de 3/4 para el metabolismo basal en mamíferos, hipótesis que posteriormente dio origen a modelos más amplios como la Metabolic Theory of Ecology (MTE), donde se plantea la existencia de principios universales capaces de explicar patrones ecológicos y fisiológicos a gran escala.
Con el tiempo; numerosos estudios posteriores comenzaron a cuestionar la universalidad de estos exponentes,es aquí donde se comienza a resaltar la importancia de tener en cuanta el efecto de los grupos taxonómicos y diferentes tipos de metabolismos. Al asumir la independencia estadística entre especies, se ignoraba que organismos emparentados evolutivamente pueden compartir características fisiológicas heredadas. Debido a esto, las estimaciones de las pendientes alométricas podrían estar sesgadas si no se incorpora explícitamente la historia filogenética en los modelos estadísticos.
De esta manera, el articulo base con el que se esta desarrollando este documento propuso reevaluar el escalamiento metabólico en mamíferos utilizando métodos comparativos filogenéticos, particularmente modelos PGLS con estimación de λ por máxima verosimilitud. La pregunta central planteada por los autores fue determinar si existe realmente un exponente metabólico universal en mamíferos o si, por el contrario, las pendientes alométricas varían entre linajes debido a la influencia de la filogenia y la historia evolutiva de las especies.
4 Materiales y métodos
- Origen de los datos
Para este documento se utilizaron datos reportados en el articulo Phylogeny and Metabolic Scaling in Mammals Capellini et al.(2010) Los dos datos clave fueron los registros de metabolismo basal (BMR) y metabolismo de campo (FMR) correspondientes a mamíferos de distinatas especies, asi como sus respectivas masas corporales.Los datos originales fueron compilados por los autores a partir de estudios fisiológicos previos y posteriormente asociados con una filogenia de mamíferos derivada del supertree desarrollado por Bininda-Emonds et al. (2007, 2008).
Para la recreacion de los datos, se usaron las mismas variables y los datos fueron transformados a escala logarítmica (log10) con el fin de poder observar de manera lineal las relaciones alométricas.
- Software utilizado
Todos los analisis fueron realizados por medio de el software R, para el manejo y visualizacion de su contenido se usaron los paquetes: tidyverse, ape, nlme, caper, ggplot2
La filogenia de mamíferos utilizada en los análisis fue procesada en formato Newick y posteriormente ajustada para coincidir con las especies presentes en las bases de datos metabólicas. Adicionalmente, se emplearon funciones de comparación filogenética para eliminar especies ausentes o inconsistencias taxonómicas entre la filogenia y los datos fisiológicos.
- Análisis estadístico (descripción)
Las relaciones alométricas entre metabolismo y masa corporal fueron evaluadas mediante modelos de regresión lineales en escala log-log. Debido a que las especies no representan observaciones estadísticamente independientes por compartir historia evolutiva, se incorporó información filogenética en los análisis mediante modelos de mínimos cuadrados generalizados filogenéticos (PGLS).
“When ML λ is intermediate between 0 and 1, therefore, both OLS and PIC are not ideal methods, because they respectively underestimate and overestimate the influence of phylogeny. Conversely, PGLS offers a more flexible approach by simultaneously estimating ML λ when testing the association between variables with regression analysis, thus accounting for the precise strength of the phylogenetic signal that the data exhibit.” (Capellini et al., 2010)
Esto es clave para entender como se va a desarrollar el análisis estadístico, en donde el parámetro λ fue estimado mediante máxima verosimilitud (Maximum Likelihood, ML), procedimiento estadístico que identifica el valor de λ que produce el mejor ajuste del modelo a los datos observados.
Para cuantificar la intensidad de la señal filogenética se estimó el parámetro λ (lambda) de Pagel mediante máxima verosimilitud. Valores de λ cercanos a 0 indican ausencia de señal filogenética y equivalen a asumir independencia entre especies, como ocurre en modelos OLS. Por el contrario, valores cercanos a 1 indican que la variación observada sigue el patrón esperado bajo evolución Browniana, produciendo resultados equivalentes a los obtenidos mediante contrastes filogenéticamente independientes (PIC).
Posteriormente, se compararon modelos con λ (PGLS), λ = 0 (OLS) y λ = 1 (PIC) mediante pruebas de razón de verosimilitud (Likelihood Ratio Test, LR), utilizando los valores de log-verosimilitud de cada modelo para identificar el ajuste estadístico más apropiado.
5 Resultados
Tabla1. Esta tabla resume la composición general de los datos que se usaron para este análisis comparativo (no los del articulo).
Así como en el articulo se trabajó con un mayor número de especies para metabolismo basal (BMR) que para metabolismo de campo (FMR), esto se debe a la mayor disponibilidad de mediciones de BMR en comparación con FMR. La mayoría de especies presentaron coincidencia con la filogenia utilizada, permitiendo incorporar adecuadamente la información evolutiva en los modelos PGLS y PIC. Asimismo, los datos incluyeron principalmente especies de Eutheria (mamíferos placentarios), mientras que Marsupialia (marsupiales) estuvo representado por un número considerablemente menor de especies. La inclusión de múltiples órdenes taxonómicos con tamaños de muestra amplios permitió realizar comparaciones alométricas entre distintos linajes de mamíferos.
Figura 1.a Presenta la relación entre el logaritmo de la tasa metabólica basal (log BMR, medida en mL O₂/h) y el logaritmo de la masa corporal (log masa, en kg) para 516 especies de mamíferos. Cada punto representa una especie, diferenciada por su subclase mediante símbolos: círculos abiertos para placentales, círculos rellenos para marsupiales y diamantes para monotremas.El parámetro b que aparece en las anotaciones es la pendiente de cada línea de regresión, es decir, el exponente alométrico que describe cuánto aumenta el metabolismo por cada unidad de aumento en la masa corporal en escala logarítmica.
Todos los grupos muestran una relación lineal positiva muy clara entre el logaritmo de la masa corporal y el logaritmo del BMR, con un R² cercano a 0.94, lo que indica que la masa corporal explica el 94% de la variación en el metabolismo basal.
Los insectívoros tienen la pendiente más baja (0.479), lo que significa que su metabolismo escala más lentamente con el tamaño que cualquier otro grupo. Esto es consistente con lo reportado en el artículo, que también encontró que los insectívoros tienen pendientes significativamente menores que todos los demás linajes. Los carnívoros y murciélagos tienen las pendientes más altas, cercanas a 0.75, mientras que roedores y placentales en general se ubican en valores intermedios alrededor de 0.69-0.70, excluyendo el exponente teórico de 3/4.
Figura 1.b Es similar a la 1.a excepto por la tasa metabólica que en este caso es de campo (log FMR, medida en kJ/día) en lugar del BMR, para 116 especies medidas en condiciones silvestres mediante el método de agua doblemente marcada. A diferencia del BMR que se mide en laboratorio bajo condiciones controladas, el FMR representa el gasto energético real diario del animal en su ambiente natural, incluyendo actividad, termorregulación y reproducción, por lo que es una medida ecológicamente más relevante.
La relación también es lineal y fuerte, con R² de 0.945, pero se nota claramente que las líneas de placentales y marsupiales se separan más que en la figura 1a, especialmente a partir de masas corporales medianas hacia arriba.
La pendiente de marsupiales (0.596) es notablemente más baja que la de placentales (0.761), lo que indica que en condiciones de vida silvestre los marsupiales gastan proporcionalmente menos energía a medida que aumentan de tamaño comparado con los placentales. Esta diferencia es biológicamente relevante porque sugiere que los marsupiales tienen una menor capacidad de incrementar su gasto energético diario con el tamaño, posiblemente relacionado con su menor tasa metabólica basal y menor capacidad de termorregulación. El artículo reporta esta misma diferencia como estadísticamente significativa (t = 3.92, p < 0.001), siendo la única comparación donde BMR y FMR difieren entre sí dentro de un mismo linaje.
Grafica 1 Se observa el contraste teniendo en cuenta que lambda representa si no hay relación en lo absoluto o hay una relación de similitud entre especies y pueden o no tratarse como independientes. Lh (Log-verosimilitud)
Tabla 2 Lh_PGLSk y Lh_otro son los valores de log-verosimilitud de cada modelo que se está comparando.LR (Likelihood Ratio) es el estadístico de prueba que cuantifica la diferencia entre los dos modelos que se comparan.
La gráfica 1 compara qué tan bien se ajusta cada modelo a los datos. A mayor log-verosimilitud (Lh), mejor ajusta el modelo.
OLS tiene una log-verosimilitud de 22.24 como se puede ver en la tabla 2, claramente inferior a los otros dos modelos. Esto confirma que ignorar completamente la filogenia produce un modelo que describe peor los datos, porque trata a especies emparentadas como si fueran independientes cuando no lo son.
PGLSk y PIC tienen exactamente el mismo valor (38.59) porque el λ recreado estimado fue 1.0 , haciendo que ambos modelos sean matemáticamente idénticos. Esto significa que en el dataset recreado la señal filogenética es máxima y corregir por filogenia mejora sustancialmente el ajuste del modelo.
6 Discusión
Los resultados obtenidos en este análisis son en gran medida consistentes con los reportados por Capellini et al. (2010), tanto en la dirección de los efectos como en las conclusiones biológicas principales.Como se puede ver las relaciones alometricas obtenidas son altamente similares con las del articulo.
(Capellini et al., 2010)
(Capellini et al., 2010)
Los análisis de regresión lineal ordinaria (OLS) replicaron con alta precisión los valores reportados en el artículo original. Para el BMR de todos los mamíferos se obtuvo una pendiente de 0.691 (IC 95%: 0.677–0.705), prácticamente idéntica al valor publicado de 0.691 (0.677–0.704). De manera similar, las pendientes por orden taxonómico mostraron el mismo patrón general descrito en el artículo: los insectívoros presentaron la pendiente más baja (0.479 vs 0.475 en el artículo), mientras que carnívoros y murciélagos mostraron los valores más altos, cercanos a 0.75. Para el FMR, la pendiente global de 0.719 fue muy cercana al valor reportado de 0.715, y la diferencia entre placentales (0.761) y marsupiales (0.596) replicó el patrón descrito por los autores, quienes reportaron esta como la única comparación significativa entre BMR y FMR dentro de un mismo linaje (t = 3.92, p < 0.001). Adicionalmente, la prueba de no linealidad mostró un incremento en R² de apenas 0.0037 al agregar el término cuadrático, resultado casi idéntico al 0.003 reportado en el artículo, confirmando que el modelo lineal es suficiente para describir la relación entre metabolismo y masa corporal en escala log-log.
Detodas formas hubo diferencias entre estos dos trabajos;mientras que Capellini et al. (2010) reportaron un λ estimado de 0.85 para el BMR de todos los mamíferos, en este análisis se obtuvo λ = 1.00, lo que produjo que PGLSk y PIC convergieran al mismo resultado (b = 0.701). Esta discrepancia tiene tres explicaciones principales.
En primer lugar, la pérdida de 64 especies al cruzar la base de datos con la filogenia de Bininda-Emonds et al. (2007) redujo el dataset de 580 a 516 especies. Las especies no coincidentes probablemente pertenecen a grupos taxonómicos específicos cuya ausencia altera la estructura filogenética del dataset y por tanto la señal filogenética que el modelo detecta.
En segundo lugar, el artículo original utilizó BayesTraits para la estimación de λ, mientras que en este análisis se empleó el paquete caper en R. Aunque ambos implementan el modelo de Pagel (1997, 1999), difieren en sus algoritmos de optimización numérica, lo que puede producir estimaciones distintas de λ especialmente cuando el valor verdadero se encuentra cerca de los límites del parámetro.
Finalmente las diferencias en nomenclatura taxonómica entre la base de datos y la filogenia, producto de actualizaciones sistemáticas desde 2007, impidieron la coincidencia de un 11% de las especies, lo que pudo haber sesgado la representación de algunos linajes en el análisis filogenético.
7 Conclusión
Aunque hubo algunas diferencias, la conclusión central del articulo se mantiene,ignorar la filogenia produce estimaciones erróneas de los exponentes alométricos, tal como argumentan Capellini et al. (2010) al señalar que los estudios previos que reportaban predominancia del exponente 3/4 cometían precisamente este error estadístico. La no significancia de la comparación PGLSk vs PIC (LR = 0.01, p = 0.938) es consistente con el λ = 1 estimado y no contradice la conclusión biológica, sino que refleja una limitación técnica de la replicación relacionada con el número de especies disponibles y el software utilizado.
En conjunto, los resultados de esta replicación apoyan las dos conclusiones principales de Capellini et al. (2010): que los exponentes alométricos varían entre linajes de mamíferos y que ninguno de los valores teóricos propuestos (0.66 o 0.75) es universalmente aplicable, y que el control estadístico por filogenia es indispensable para obtener estimaciones correctas de dichos exponentes. Las diferencias numéricas observadas entre este análisis y el artículo original son metodológicamente explicables y no alteran la interpretación biológica de los resultados, lo que valida la robustez de las conclusiones del estudio original.
8 Código reproducible
En esta sección esta el código usado para las imágenes que se descargaron y agregaron al articulo, así como algunos datos que aunque puede que no se usaron explicitamente en el análisis, fueron claves para la obtención de los datos que si fueron discutidos en el documento. Para la tabla 2 y gráfica 1
(
LR_A_logverosimilitud.pngLR_B_tabla.png)
Estas fueron descargadas directamente como imágenes a mi computador por lo que no grafican en esta sección pero se ven en el documento.
# =============================================================================
# REPLICACIÓN: Capellini et al. (2010)
# "Phylogeny and metabolic scaling in mammals"
# Ecology, 91(9): 2783-2793
# =============================================================================
# Paquetes necesarios:
# install.packages(c("tidyverse", "ape", "nlme", "ggplot2", "caper"))
# El paquete 'caper' incluye funciones PGLS similares a BayesTraits
library(tidyverse)
library(ape)
library(nlme)
library(caper)
library(ggplot2)
# -----------------------------------------------------------------------------
# 1. CARGA Y LIMPIEZA DE DATOS
# -----------------------------------------------------------------------------
datos_raw <- read.table(
"C:/Users/mcami/OneDrive/Desktop/data.txt", # ajusta la ruta si es necesario
header = TRUE,
sep = "\t",
stringsAsFactors = FALSE,
check.names = FALSE
)
# Renombrar columnas para facilitar el manejo
colnames(datos_raw) <- c("Subclass", "Order", "Species",
"BMR", "Mass_BMR",
"Ref_BMR",
"FMR", "Mass_FMR",
"Ref_FMR")
# Reemplazar -9999 por NA (código de datos faltantes en el archivo)
datos_raw[datos_raw == -9999] <- NA
# ---- Dataset BMR ----
# Incluye solo especies con BMR válido
# Masa en gramos → convertir a kg para log-transformación coherente con el artículo
bmr_data <- datos_raw %>%
filter(!is.na(BMR), !is.na(Mass_BMR)) %>%
mutate(
log_BMR = log10(BMR), # log10(mL O2/h)
log_Mass_BMR = log10(Mass_BMR / 1000) # log10(kg)
)
cat("N especies con BMR:", nrow(bmr_data), "\n") # esperado ~580
# ---- Dataset FMR ----
fmr_data <- datos_raw %>%
filter(!is.na(FMR), !is.na(Mass_FMR)) %>%
mutate(
log_FMR = log10(FMR), # log10(kJ/d)
log_Mass_FMR = log10(Mass_FMR / 1000) # log10(kg)
)
cat("N especies con FMR:", nrow(fmr_data), "\n") # esperado ~119
# -----------------------------------------------------------------------------
# 2. REGRESIÓN OLS (k = 0, equivalente a species-level analysis)
# -----------------------------------------------------------------------------
# Esta es la regresión lineal ordinaria que el artículo critica por ignorar
# la filogenia. Se usa como punto de comparación (Tabla 5 del artículo).
# -- OLS para BMR (todos los mamíferos) --
ols_bmr_all <- lm(log_BMR ~ log_Mass_BMR, data = bmr_data)
summary(ols_bmr_all)
# Extraer pendiente e IC 95%
coef_bmr_ols <- coef(ols_bmr_all)["log_Mass_BMR"]
ci_bmr_ols <- confint(ols_bmr_all)["log_Mass_BMR", ]
cat("\n--- OLS BMR todos los mamíferos ---\n")
cat("Pendiente (b):", round(coef_bmr_ols, 3), "\n")
cat("IC 95%:", round(ci_bmr_ols, 3), "\n")
# Artículo reporta: b = 0.691 (0.677-0.704)
# -- OLS para FMR (todos los mamíferos) --
ols_fmr_all <- lm(log_FMR ~ log_Mass_FMR, data = fmr_data)
summary(ols_fmr_all)
coef_fmr_ols <- coef(ols_fmr_all)["log_Mass_FMR"]
ci_fmr_ols <- confint(ols_fmr_all)["log_Mass_FMR", ]
cat("\n--- OLS FMR todos los mamíferos ---\n")
cat("Pendiente (b):", round(coef_fmr_ols, 3), "\n")
cat("IC 95%:", round(ci_fmr_ols, 3), "\n")
# Artículo reporta: b = 0.715 (0.684-0.745)
# -- OLS por subclase --
subclases_bmr <- c("Eutheria", "Marsupialia")
for (sc in subclases_bmr) {
d <- bmr_data %>% filter(Subclass == sc)
m <- lm(log_BMR ~ log_Mass_BMR, data = d)
b <- round(coef(m)["log_Mass_BMR"], 3)
ic <- round(confint(m)["log_Mass_BMR", ], 3)
cat("\nOLS BMR", sc, "| b =", b, "| IC:", ic, "| n =", nrow(d), "\n")
}
# -- OLS por orden (órdenes con n > 30) --
ordenes_grandes_bmr <- bmr_data %>%
count(Order) %>%
filter(n >= 30) %>%
pull(Order)
cat("\nÓrdenes con n >= 30 para BMR:", ordenes_grandes_bmr, "\n")
for (ord in ordenes_grandes_bmr) {
d <- bmr_data %>% filter(Order == ord)
m <- lm(log_BMR ~ log_Mass_BMR, data = d)
b <- round(coef(m)["log_Mass_BMR"], 3)
ic <- round(confint(m)["log_Mass_BMR", ], 3)
cat("OLS BMR", ord, "| b =", b, "| IC:", ic, "| n =", nrow(d), "\n")
}
# -----------------------------------------------------------------------------
# 3. PGLS CON ESTIMACIÓN DE LAMBDA (equivalente a PGLSk del artículo)
# -----------------------------------------------------------------------------
# NOTA IMPORTANTE:
# El artículo original usa BayesTraits con la filogenia de Bininda-Emonds et al.
# (2007, 2008). Aquí se muestra la estructura del código PGLS usando el paquete
# 'caper', que implementa el mismo modelo estadístico (Freckleton et al. 2002).
# Para una réplica exacta se necesita cargar el supertree de mamíferos.
# ---- 3a. Cargar filogenia (requiere archivo externo) ----
# La filogenia de Bininda-Emonds está disponible en:
# https://doi.org/10.1038/nature05634
# Descargar el archivo .nex o .tre y cargarlo así:
#
# filogenia <- read.nexus("bininda_emonds_2007.nex")
# -- O en formato Newick:
# filogenia <- read.tree("bininda_emonds_2007.tre")
#
# Si no se tiene la filogenia, se puede usar ape::rtree() para demostración:
#install.packages("PhyloOrchard", repos = "http://R-Forge.R-project.org")
library(PhyloOrchard)
data(BinindaEmondsEtAl2007)
filogenia_completa <- BinindaEmondsEtAl2007[[1]] # best estimate
# Estandarizar nombres: espacios → guiones bajos
bmr_data$Species <- gsub(" ", "_", bmr_data$Species)
fmr_data$Species <- gsub(" ", "_", fmr_data$Species)
# Podar el árbol a las especies de tu dataset
especies_bmr <- bmr_data$Species
coinciden_bmr <- especies_bmr[especies_bmr %in% filogenia_completa$tip.label]
cat("Especies BMR en el árbol:", length(coinciden_bmr),
"de", length(especies_bmr), "\n")
filogenia_demo <- keep.tip(filogenia_completa, coinciden_bmr)
bmr_data <- bmr_data %>% filter(Species %in% coinciden_bmr)
cat("Especies BMR que coinciden con el árbol:", length(coinciden_bmr), "\n")
# ---- 3b. Crear objeto comparative.data (requerido por caper) ----
# Combina los datos con la filogenia
comp_data_bmr <- comparative.data(
phy = filogenia_demo,
data = bmr_data,
names.col = "Species",
vcv = TRUE,
warn.dropped = FALSE
)
# ---- 3c. PGLS con lambda estimado por máxima verosimilitud ----
# lambda = "ML" → estima λ simultáneamente con la pendiente (= PGLSk del artículo)
# lambda = 0 → equivale a OLS
# lambda = 1 → equivale a PIC
pgls_bmr_lambda_ML <- pgls(
log_BMR ~ log_Mass_BMR,
data = comp_data_bmr,
lambda = "ML" # ← estima λ por máxima verosimilitud
)
summary(pgls_bmr_lambda_ML)
# Lambda estimado
lambda_ML_bmr <- pgls_bmr_lambda_ML$param["lambda"]
cat("\nLambda ML para BMR (todos):", round(lambda_ML_bmr, 3), "\n")
# Artículo reporta: λ = 0.85
# Pendiente e IC 95%
b_pgls_bmr <- coef(pgls_bmr_lambda_ML)["log_Mass_BMR"]
# Extraer IC manualmente porque confint() no funciona bien con pgls
se_pgls_bmr <- summary(pgls_bmr_lambda_ML)$coefficients["log_Mass_BMR", "Std. Error"]
ic_pgls_bmr <- c(
b_pgls_bmr - 1.96 * se_pgls_bmr,
b_pgls_bmr + 1.96 * se_pgls_bmr
)
names(ic_pgls_bmr) <- c("2.5%", "97.5%")
cat("Pendiente PGLSk BMR:", round(b_pgls_bmr, 3), "\n")
cat("IC 95%:", round(ic_pgls_bmr, 3), "\n")
# Artículo reporta: b = 0.718 (0.697-0.738)
# ---- 3d. Comparar modelos con Likelihood Ratio Test ----
# PGLSk vs OLS (λ = 0)
pgls_bmr_lambda0 <- pgls(log_BMR ~ log_Mass_BMR, data = comp_data_bmr, lambda = 0.0001)
pgls_bmr_lambda1 <- pgls(log_BMR ~ log_Mass_BMR, data = comp_data_bmr, lambda = 0.9999)
# LR test: LR = 2 * (Lh_mejor - Lh_peor), distribuido como chi² con 1 gl
LR_vs_OLS <- 2 * (logLik(pgls_bmr_lambda_ML) - logLik(pgls_bmr_lambda0))
LR_vs_PIC <- 2 * (logLik(pgls_bmr_lambda_ML) - logLik(pgls_bmr_lambda1))
p_vs_OLS <- pchisq(as.numeric(LR_vs_OLS), df = 1, lower.tail = FALSE)
p_vs_PIC <- pchisq(as.numeric(LR_vs_PIC), df = 1, lower.tail = FALSE)
cat("\n--- LR Test BMR: PGLSk vs OLS ---\n")
cat("LR =", round(as.numeric(LR_vs_OLS), 2), "| p =", round(p_vs_OLS, 5), "\n")
# p < 0.05 → PGLSk ajusta significativamente mejor que OLS
cat("--- LR Test BMR: PGLSk vs PIC ---\n")
cat("LR =", round(as.numeric(LR_vs_PIC), 2), "| p =", round(p_vs_PIC, 5), "\n")
# -----------------------------------------------------------------------------
# LR TEST - Tres figuras separadas y limpias
# -----------------------------------------------------------------------------
library(ggplot2)
library(gridExtra)
library(grid)
# Datos base
lh_pgls <- as.numeric(logLik(pgls_bmr_lambda_ML))
lh_ols <- as.numeric(logLik(pgls_bmr_lambda0))
lh_pic <- as.numeric(logLik(pgls_bmr_lambda1))
LR_pgls_vs_ols <- 2 * (lh_pgls - lh_ols)
LR_pgls_vs_pic <- 2 * (lh_pgls - lh_pic)
p_pgls_vs_ols <- pchisq(LR_pgls_vs_ols, df = 1, lower.tail = FALSE)
p_pgls_vs_pic <- pchisq(LR_pgls_vs_pic, df = 1, lower.tail = FALSE)
valor_critico <- qchisq(0.95, df = 1)
# =============================================================
# FIGURA A: Log-verosimilitud
# =============================================================
datos_lh <- data.frame(
Modelo = factor(c("OLS (λ=0)", "PGLSk (λ=ML)", "PIC (λ=1)"),
levels = c("OLS (λ=0)", "PGLSk (λ=ML)", "PIC (λ=1)")),
Lh = c(lh_ols, lh_pgls, lh_pic),
Es_mejor = c(FALSE, TRUE, FALSE)
)
figA <- ggplot(datos_lh, aes(x = Modelo, y = Lh, fill = Es_mejor)) +
geom_bar(stat = "identity", width = 0.5,
color = "black", linewidth = 0.5) +
geom_text(aes(label = round(Lh, 2)),
vjust = -0.8, size = 4, fontface = "bold") +
scale_fill_manual(values = c("FALSE" = "gray80",
"TRUE" = "#27ae60")) +
scale_y_continuous(expand = expansion(mult = c(0, 0.2))) +
labs(
title = "Log-verosimilitud por modelo",
x = "Modelo",
y = "Log-verosimilitud (Lh)"
) +
theme_classic(base_size = 13) +
theme(
legend.position = "none",
plot.title = element_text(face = "bold", hjust = 0.5, size = 13),
axis.text = element_text(size = 11)
)
ggsave("LR_A_logverosimilitud.png",
figA, width = 6, height = 5, dpi = 300)
# =============================================================
# FIGURA B: Tabla LR test
# =============================================================
tabla_df <- data.frame(
Comparacion = c("PGLSk vs OLS", "PGLSk vs PIC"),
Lh_PGLSk = c(round(lh_pgls, 2), round(lh_pgls, 2)),
Lh_otro = c(round(lh_ols, 2), round(lh_pic, 2)),
LR = c(round(LR_pgls_vs_ols, 2), round(LR_pgls_vs_pic, 2)),
P_valor = c("<0.0001", round(p_pgls_vs_pic, 4)),
Conclusion = c("PGLSk superior", "Sin diferencia")
)
colores_filas <- c("#d5f5e3", "#fadbd8")
figB <- tableGrob(
tabla_df,
rows = NULL,
theme = ttheme_default(
colhead = list(
fg_params = list(fontface = "bold",
fontsize = 11,
col = "white"),
bg_params = list(fill = "gray20")
),
core = list(
fg_params = list(fontsize = 10),
bg_params = list(
fill = matrix(
c(colores_filas[1], colores_filas[2]),
nrow = 2, ncol = ncol(tabla_df)
)
)
),
padding = unit(c(8, 5), "mm")
)
)
png("LR_B_tabla.png", width = 900, height = 200, res = 120)
grid.newpage()
grid.draw(figB)
dev.off()
# =============================================================
# FIGURA C: Barras LR
# =============================================================
datos_lr <- data.frame(
Comparacion = factor(
c("PGLSk vs OLS", "PGLSk vs PIC"),
levels = c("PGLSk vs OLS", "PGLSk vs PIC")
),
LR = c(LR_pgls_vs_ols, LR_pgls_vs_pic),
Significativo = c(TRUE, FALSE),
P_label = c("p < 0.0001",
paste0("p = ", round(p_pgls_vs_pic, 4)))
)
figC <- ggplot(datos_lr,
aes(x = Comparacion, y = LR,
fill = Significativo)) +
geom_bar(stat = "identity", width = 0.4,
color = "black", linewidth = 0.5) +
geom_text(aes(label = paste0("LR = ", round(LR, 2))),
vjust = -0.5, size = 4, fontface = "bold") +
geom_text(aes(label = P_label),
vjust = -2.0, size = 3.5, color = "gray30") +
geom_hline(yintercept = valor_critico,
linetype = "dashed",
color = "red", linewidth = 0.9) +
annotate("text",
x = 1.5,
y = valor_critico + 1.5,
label = paste0("Valor crítico χ²(1) = ",
round(valor_critico, 2)),
color = "red", size = 3.5, hjust = 0.5) +
scale_fill_manual(
values = c("FALSE" = "gray80", "TRUE" = "#2980b9"),
labels = c("No significativo", "Significativo")
) +
scale_y_continuous(expand = expansion(mult = c(0, 0.3))) +
labs(
title = "Estadístico LR por comparación",
x = "Comparación",
y = "Estadístico LR",
fill = ""
) +
theme_classic(base_size = 13) +
theme(
legend.position = "bottom",
plot.title = element_text(face = "bold",
hjust = 0.5, size = 13),
axis.text = element_text(size = 11)
)
ggsave("LR_C_estadistico.png",
figC, width = 6, height = 5, dpi = 300)
cat("Figuras guardadas:\n")
cat(" - LR_A_logverosimilitud.png\n")
cat(" - LR_B_tabla.png\n")
cat(" - LR_C_estadistico.png\n")
# ---- 3e. PIC (λ = 1, Phylogenetically Independent Contrasts) ----
pgls_bmr_PIC <- pgls(log_BMR ~ log_Mass_BMR, data = comp_data_bmr, lambda = 1)
b_pic_bmr <- coef(pgls_bmr_PIC)["log_Mass_BMR"]
se_pic_bmr <- summary(pgls_bmr_PIC)$coefficients["log_Mass_BMR", "Std. Error"]
ic_pic_bmr <- c(
b_pic_bmr - 1.96 * se_pic_bmr,
b_pic_bmr + 1.96 * se_pic_bmr
)
names(ic_pic_bmr) <- c("2.5%", "97.5%")
cat("\nPIC BMR (λ=1) | b =", round(b_pic_bmr, 3),
"| IC:", round(ic_pic_bmr, 3), "\n")
# Artículo reporta: b = 0.743 (0.717-0.769)
# ---- 3f. Resumen de las tres pendientes (Tabla 5 del artículo) ----
tabla_comparacion <- data.frame(
Metodo = c("OLS (k=0)", "PGLSk (k=ML)", "PIC (k=1)"),
Pendiente = c(round(coef_bmr_ols, 3),
round(b_pgls_bmr, 3),
round(b_pic_bmr, 3)),
IC_inf = c(round(ci_bmr_ols[1], 3),
round(ic_pgls_bmr[1], 3),
round(ic_pic_bmr[1], 3)),
IC_sup = c(round(ci_bmr_ols[2], 3),
round(ic_pgls_bmr[2], 3),
round(ic_pic_bmr[2], 3))
)
print(tabla_comparacion)
# -----------------------------------------------------------------------------
# 6. FIGURAS PUBLICABLES - Réplica de Fig. 1a y 1b del artículo
# -----------------------------------------------------------------------------
# Valores de pendientes por orden para anotar en la figura (del artículo)
b_rodentia <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Order == "Rodentia")))["log_Mass_BMR"], 3)
b_chiroptera <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Order == "Chiroptera")))["log_Mass_BMR"], 3)
b_carnivora <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Order == "Carnivora")))["log_Mass_BMR"], 3)
b_lypothyphla <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Order == "Lypothyphla")))["log_Mass_BMR"], 3)
b_eutheria <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Subclass == "Eutheria")))["log_Mass_BMR"], 3)
b_marsupialia <- round(coef(lm(log_BMR ~ log_Mass_BMR, bmr_data %>% filter(Subclass == "Marsupialia")))["log_Mass_BMR"], 3)
# Calcular pendientes BMR antes de la figura
b_rodentia <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Order == "Rodentia")))["log_Mass_BMR"], 3)
b_chiroptera <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Order == "Chiroptera")))["log_Mass_BMR"], 3)
b_carnivora <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Order == "Carnivora")))["log_Mass_BMR"], 3)
b_lypothyphla <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Order == "Lypothyphla")))["log_Mass_BMR"], 3)
b_eutheria <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Subclass == "Eutheria")))["log_Mass_BMR"], 3)
b_marsupialia <- round(coef(lm(log_BMR ~ log_Mass_BMR,
bmr_data %>% filter(Subclass == "Marsupialia")))["log_Mass_BMR"], 3)
# ---- Figura 1a: BMR vs masa corporal ----
fig1a <- ggplot(bmr_data, aes(x = log_Mass_BMR, y = log_BMR)) +
geom_point(aes(shape = Subclass),
color = "gray30", alpha = 0.6, size = 1.5) +
scale_shape_manual(
values = c("Eutheria" = 1, "Marsupialia" = 16, "Monotremata" = 5),
labels = c("Placentales", "Marsupiales", "Monotremas")
) +
# Líneas por orden
geom_smooth(data = bmr_data %>% filter(Order == "Rodentia"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "#E69F00", linewidth = 0.9) +
geom_smooth(data = bmr_data %>% filter(Order == "Chiroptera"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "#0072B2", linewidth = 0.9) +
geom_smooth(data = bmr_data %>% filter(Order == "Carnivora"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "#009E73", linewidth = 0.9) +
geom_smooth(data = bmr_data %>% filter(Order == "Lypothyphla"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "#D55E00", linewidth = 0.9) +
# Líneas por subclase
geom_smooth(data = bmr_data %>% filter(Subclass == "Eutheria"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "black", linewidth = 0.8, linetype = "solid") +
geom_smooth(data = bmr_data %>% filter(Subclass == "Marsupialia"),
aes(x = log_Mass_BMR, y = log_BMR),
method = "lm", se = FALSE,
color = "black", linewidth = 0.8, linetype = "dashed") +
# Etiquetas de órdenes en la DERECHA del gráfico, bien separadas
annotate("text", x = 2.5, y = 2.8,
label = paste0("Rodents: b = ", b_rodentia),
color = "#E69F00", hjust = 1, size = 3.0) +
annotate("text", x = 2.5, y = 2.5,
label = paste0("Bats: b = ", b_chiroptera),
color = "#0072B2", hjust = 1, size = 3.0) +
annotate("text", x = 2.5, y = 2.2,
label = paste0("Carnivores: b = ", b_carnivora),
color = "#009E73", hjust = 1, size = 3.0) +
annotate("text", x = 2.5, y = 1.9,
label = paste0("Insectivores: b = ", b_lypothyphla),
color = "#D55E00", hjust = 1, size = 3.0) +
# Texto de subclases abajo a la derecha, cada línea separada
annotate("text", x = 2.5, y = 1.3,
label = "All mammals (not shown): b = 0.718",
hjust = 1, size = 2.7, color = "black") +
annotate("text", x = 2.5, y = 1.1,
label = paste0("Placentals (solid black): b = ", b_eutheria),
hjust = 1, size = 2.7, color = "black") +
annotate("text", x = 2.5, y = 0.9,
label = paste0("Marsupials (dashed black): b = ", b_marsupialia),
hjust = 1, size = 2.7, color = "black") +
# Etiqueta "a"
annotate("text", x = -2.5, y = 5.0,
label = "a", size = 5, fontface = "bold") +
labs(
x = "log(body mass)",
y = "log(BMR)",
shape = NULL
) +
# Leyenda abajo para que no tape nada
theme_classic(base_size = 12) +
theme(
legend.position = "bottom",
legend.background = element_blank(),
legend.key.size = unit(0.4, "cm"),
legend.text = element_text(size = 9),
axis.title = element_text(size = 11)
)
print(fig1a)
ggsave("figura1a_BMR.png", fig1a, width = 7, height = 6, dpi = 600)
# Calcular pendientes FMR antes de la figura
b_eutheria_fmr <- round(coef(lm(log_FMR ~ log_Mass_FMR,
fmr_data %>% filter(Subclass == "Eutheria")))["log_Mass_FMR"], 3)
b_marsupialia_fmr <- round(coef(lm(log_FMR ~ log_Mass_FMR,
fmr_data %>% filter(Subclass == "Marsupialia")))["log_Mass_FMR"], 3)
b_todos_fmr <- round(coef(lm(log_FMR ~ log_Mass_FMR,
fmr_data))["log_Mass_FMR"], 3)
# ---- Figura 1b: FMR vs masa corporal ----
fig1b <- ggplot(fmr_data, aes(x = log_Mass_FMR, y = log_FMR)) +
geom_point(aes(shape = Subclass),
color = "gray30", alpha = 0.6, size = 1.5) +
scale_shape_manual(
values = c("Eutheria" = 1, "Marsupialia" = 16, "Monotremata" = 5),
labels = c("Placentales", "Marsupiales", "Monotremas")
) +
geom_smooth(data = fmr_data %>% filter(Subclass == "Eutheria"),
aes(x = log_Mass_FMR, y = log_FMR),
method = "lm", se = FALSE,
color = "black", linewidth = 0.8, linetype = "solid") +
geom_smooth(data = fmr_data %>% filter(Subclass == "Marsupialia"),
aes(x = log_Mass_FMR, y = log_FMR),
method = "lm", se = FALSE,
color = "black", linewidth = 0.8, linetype = "dashed") +
# Texto bien separado línea por línea
annotate("text", x = 2.0, y = 1.9,
label = "All mammals (not shown): b = 0.697",
hjust = 1, size = 2.7, color = "black") +
annotate("text", x = 2.0, y = 1.7,
label = paste0("Placentals (solid black): b = ", b_eutheria_fmr),
hjust = 1, size = 2.7, color = "black") +
annotate("text", x = 2.0, y = 1.5,
label = paste0("Marsupials (dashed black): b = ", b_marsupialia_fmr),
hjust = 1, size = 2.7, color = "black") +
annotate("text", x = -2.2, y = 4.7,
label = "b", size = 5, fontface = "bold") +
labs(
x = "log(body mass)",
y = "log(FMR)",
shape = NULL
) +
theme_classic(base_size = 12) +
theme(
legend.position = "bottom",
legend.background = element_blank(),
legend.key.size = unit(0.4, "cm"),
legend.text = element_text(size = 9),
axis.title = element_text(size = 11)
)
print(fig1b)
ggsave("figura1b_FMR.png", fig1b, width = 7, height = 6, dpi = 600)
# ============================================================
# RESUMEN GENERAL DEL DATASET
# ============================================================
library(dplyr)
library(knitr)
library(kableExtra)
tabla_resumen <- data.frame(
Categoria = c(
"Especies con datos de BMR",
"Especies con datos de FMR",
"Especies en filogenia completa",
"Especies BMR coincidentes con filogenia",
"Especies BMR excluidas",
"Eutheria (BMR)",
"Marsupialia (BMR)",
"Órdenes BMR con n ≥ 30"
),
N = c(
nrow(bmr_data),
nrow(fmr_data),
length(filogenia_completa$tip.label),
length(coinciden_bmr),
nrow(bmr_data) - length(coinciden_bmr),
nrow(bmr_data %>% filter(Subclass == "Eutheria")),
nrow(bmr_data %>% filter(Subclass == "Marsupialia")),
length(ordenes_grandes_bmr)
)
)
# Mostrar tabla
kable(
tabla_resumen,
caption = "Tabla 1. Resumen general de los datos utilizados en los análisis."
) %>%
kable_styling(
full_width = FALSE,
bootstrap_options = c("striped", "hover", "condensed")
)9 Referencias
9.1 Referencias de software y paquetes utilizados
R Core Team (2025). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.
9.1.1 tidyverse
Wickham H, Averick M, Bryan J, Chang W, McGowan LD, François R, Grolemund G, Hayes A, Henry L, Hester J, Kuhn M, Pedersen TL, Miller E, Bache SM, Müller K, Ooms J, Robinson D, Seidel DP, Spinu V, Takahashi K, Vaughan D, Wilke C, Woo K, Yutani H (2019). “Welcome to the tidyverse.” Journal of Open Source Software, 4(43), 1686. doi:10.21105/joss.01686 https://doi.org/10.21105/joss.01686.
9.1.2 ape
Paradis E, Schliep K (2019). “ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R.” Bioinformatics, 35, 526-528. doi:10.1093/bioinformatics/bty633 https://doi.org/10.1093/bioinformatics/bty633.
9.1.3 nlme
Pinheiro J, Bates D, R Core Team (2026). nlme: Linear and Nonlinear Mixed Effects Models. doi:10.32614/CRAN.package.nlme https://doi.org/10.32614/CRAN.package.nlme, R package version 3.1-169, https://CRAN.R-project.org/package=nlme. Pinheiro JC, Bates DM (2000). Mixed-Effects Models in S and S-PLUS. Springer, New York. doi:10.1007/b98882 https://doi.org/10.1007/b98882.
9.1.4 caper
Orme D, Freckleton R, Thomas G, Petzoldt T, Fritz S, Isaac N, Pearse W (2025). caper: Comparative Analyses of Phylogenetics and Evolution in R. doi:10.32614/CRAN.package.caper https://doi.org/10.32614/CRAN.package.caper, R package version 1.0.4, https://CRAN.R-project.org/package=caper.
9.1.5 ggplot2
Wickham H (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York. ISBN 978-3-319-24277-4, https://ggplot2.tidyverse.org.
9.1.6 PhyloOrchard
O’Meara BC, Harmon L, Eastman J (2013). PhyloOrchard: Important and/or useful phylogenetic datasets. R package version 0.3/r66, https://R-Forge.R-project.org/projects/phyloorchard/.
9.1.7 knitr
Xie Y (2025). knitr: A General-Purpose Package for Dynamic Report Generation in R. R package version 1.51, https://yihui.org/knitr/. Xie Y (2015). Dynamic Documents with R and knitr, 2nd edition. Chapman and Hall/CRC, Boca Raton, Florida. ISBN 978-1498716963, https://yihui.org/knitr/. Xie Y (2014). “knitr: A Comprehensive Tool for Reproducible Research in R.” In Stodden V, Leisch F, Peng RD (eds.), Implementing Reproducible Computational Research. Chapman and Hall/CRC. ISBN 978-1466561595.
9.1.8 kableExtra
Zhu H (2024). kableExtra: Construct Complex Table with ‘kable’ and Pipe Syntax. doi:10.32614/CRAN.package.kableExtra https://doi.org/10.32614/CRAN.package.kableExtra, R package version 1.4.0, https://CRAN.R-project.org/package=kableExtra.
10 Referencia del articulo
Capellini, I., Venditti, C., & Barton, R. A. (2010). Phylogeny and metabolic scaling in mammals. Ecology, 91(9), 2783–2793. https://doi.org/10.1890/09-0817.1