Variables LMCD y C12CD
América Georgette Villalobos Martínez
2025-11-26
R Markdown
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
summary(cars)## speed dist
## Min. : 4.0 Min. : 2.00
## 1st Qu.:12.0 1st Qu.: 26.00
## Median :15.0 Median : 36.00
## Mean :15.4 Mean : 42.98
## 3rd Qu.:19.0 3rd Qu.: 56.00
## Max. :25.0 Max. :120.00Including Plots
You can also embed plots, for example:
Note that the echo = FALSE parameter was added to the
code chunk to prevent printing of the R code that generated the
plot.
Defino mi directorio de trabajo
setwd(“C:/Users/ameri/OneDrive/Desktop/ESTADISTICA FORENSE”) ##Abriendo paquete pacman library(pacman) ## El archivo esta en formato SPSS, lo abrimos mediante la libreria haven ## MASS permite analisis discriminante p_load(haven,dplyr,ggplot2,MASS) p_load(haven,dplyr,ggplot2,tinytex,tidyr, GGally) Hombro <- read_sav(“Datos hombro.sav”)
Definimos como factor la variable sexoN
Hombro\(sexoN <- factor(Hombro\)sexoN, levels = c(1, 2), labels = c(“Hombre”, “Mujer”)) ##Frecuencias de sexo table (Hombro$sexoN)
Omitimos los NA de la tabla
HombroS <- na.omit(Hombro[, c(“sexoN”, “LMCD”,“C12CD”)])
Para resumir por sexo la información de Longitud Maxima de clavicula derecha
res_lmcd <- Hombro %>%
group_by(sexoN) %>%
summarise(
n = sum(!is.na(LMCD)),
media = mean(LMCD, na.rm = TRUE),
sd = sd(LMCD, na.rm = TRUE)
) %>%
mutate(across(c(media, sd), ~round(.x, 2)))
res_lmcd
Para resumir por sexo la información de Circunferencia media de clavicula derecha
res_C12CD <- HombroS %>%
group_by(sexoN) %>%
summarise(
n = sum(!is.na(C12CD)),
media = mean(C12CD, na.rm = TRUE),
sd = sd(C12CD, na.rm = TRUE)
) %>%
mutate(across(c(media, sd), ~round(.x, 2)))
res_C12CD
## Gráfica comparativa ggplot(HombroS, aes(x = LMCD, fill = sexoN)) + geom_density(alpha = 0.5) + labs( title = “Gráfica 1. Longitud maxima clavicula derecha por sexo”, x = “Longitud clavicula derecha por sexo (mm)”, y = “Densidad” ) + theme_minimal()
## Gráfica comparativa ggplot(HombroS, aes(x = LMCD, fill = sexoN)) + geom_density(alpha = 0.5) + labs( title = “Gráfica 1. Longitud maxima clavicula derecha por sexo”, x = “Longitud clavicula derecha por sexo (mm)”, y = “Densidad” ) + theme_minimal() ggplot(HombroS, aes(x = C12CD, fill = sexoN)) + geom_density(alpha = 0.5) + labs( title = “Gráfica 2. Circunferencia en punto medio de la clavicula derecha”, x = “Circunferencia en punto medio de la clavicula derecha (mm)”, y = “Densidad” ) + theme_minimal()
1. Comparaciones por sexo entre las variables
## Realizaremos la prueba t para comparar las medias
t.test(LMCD ~ sexoN, data = HombroS, var.equal = TRUE)
t.test(LMCD ~ sexoN, data = HombroS, var.equal = FALSE)
t.test(C12CD ~ sexoN, data = HombroS, var.equal = TRUE)
t.test(C12CD ~ sexoN, data = HombroS, var.equal = FALSE)
los valores son parametricos por lo que no es necesario hacer la prueba U
2. Construir una función discriminante para LMCD
Ajustar Discriminante Lineal
# Eliminar NAs de forma explícita
HombroS <- na.omit(Hombro[, c("sexoN", "LMCD")])
# Ajustar el modelo con la base depurada
## El modelo discriminante es D=a*LMCD+b
lda1 <- lda(sexoN ~ LMCD, data = HombroS)
a <- coef(lda1) ## El factor que multiplica a LMCD
pred0 <- predict(lda1,
newdata = data.frame(LMCD = 0))
b <- pred0$x # valor de D cuando LMCD=0 es b
pred <- predict(lda1) ## Predicción de sexo para todos los valores de
# medias de la función discriminante por grupo
centroide_H <- mean(pred$x[HombroS$sexoN == "Hombre"])
centroide_M <- mean(pred$x[HombroS$sexoN == "Mujer"])
# Punto de corte (promedio de centroides, priors iguales)
cutoff <- mean(c(centroide_H, centroide_M))
## El modelo discriminante es
cat("Función discriminante: D(x) = ", round(a, 4), " * LMCD + ", round(b, 4), "\n")
Función discriminante: D(x) = 0.1089 * LMCD + -15.5587
cat(“Punto de corte=”, round(cutoff, 4), “Si D>”,round(cutoff, 4), “es Hombre”)
Función discriminante: D(x) = 0.1089 * LMCD + -15.5587
Punto de corte= -0.1791 Si D> -0.1791 es Hombre
# Tabla cruzada
tabla_clas <- table(Observado = HombroS$sexoN,
Predicho = pred$class)
tabla_clas
prop_clas <- sum(diag(tabla_clas)) / sum(tabla_clas)*100 cat(“El porcentaje de clasificación correcta es”,round(prop_clas,1),“%”)
El porcentaje de clasificación correcta es 90 %
Ejemplo de aplicación
Valores del modelo obtenidos
a <- 0.1089 # coeficiente de LMCD
b <- -15.5587 # intercepto
cutoff <- -0.1791 # punto de corte
Valor dado
LMCD_nuevo <- 140
Calcular la función discriminante
D_x <- a * LMCD_nuevo + b
Mostrar el valor calculado
cat("El valor discriminante para LMCD=140 es:", round(D_x, 4), "\n")
##El valor discriminante para LMCD=140 es: -0.3059
Clasificación según el punto de corte
if (D_x > cutoff) {
cat("Clasificación: Hombre\n")
} else {
cat("Clasificación: Mujer\n")
}Se clasifica como mujer cuando LMCD= 130
3. Construir un modelo discriminante para las dos varables
##Ajustar Disciminante Lineal para dos variables # Eliminar NAs de forma explícita Hombro_sinNA <- na.omit(Hombro[, c(“sexoN”, “LMCD”,“C12CD”)]) # Ajustar el modelo con la base depurada ## El modelo discriminante que queremos ajustar es D=aLMCD+bC12CD+c lda2 <- lda(sexoN ~ LMCD+C12CD, data = Hombro_sinNA) ##Ajustar Disciminante Lineal # Obtener los coeficientes a y b a1 <- coef(lda2)[1] a2 <- coef(lda2)[2] lda2
Predicción de valores mediante modelo
pred2 <- predict(lda2)
Matriz de confusión
table(Real = HombroS$sexoN, Predicho = pred2$class)
Porcentaje de acierto
mean(pred2$class == HombroS$sexoN)*100
valor de D cuando LMCD=0 y C12CD=0 (intercepto)
b <- predict(lda2,
newdata = data.frame(LMCD = 0, C12CD = 0))$x
b
centroide_H <- mean(pred2\(x[Hombro_sinNA\)sexoN == “Hombre”]) centroide_M <- mean(pred2\(x[Hombro_sinNA\)sexoN == “Mujer”]) cutoff <- mean(c(centroide_H, centroide_M)) cutoff
cat(“Función discriminante: D(x) =”, round(a1, 4), ” * LMHD + “, round(a2, 4),”*ABHD +“,round(b, 4),”“)
cat(“Punto de corte=”, round(cutoff, 4), “Si D>”,round(cutoff, 4), “es Hombre”)
Función discriminante: D(x) = -0.0727 * LMHD + -0.256 *ABHD + 18.9578
Punto de corte= 0.236 Si D> 0.236 es Hombre
ggplot(Hombro_sinNA, aes(x = LMCD, y = C12CD, color = sexoN)) +
geom_point(size = 3) +
theme_minimal() +
labs(title = "Discriminante lineal: LMCDD vs C12CD")
a1 <- lda2\(scaling["LMCD", 1] a2 <- lda2\)scaling[“C12CD”, 1] b <- predict(lda2, newdata = data.frame(LMCD = 0, C12CD = 0))$x - cutoff
Ecuación de la recta: a1LMCD + a2C12CD + b = 0
=> C12CD = -(a1/a2)*LMCD - b/a2
a1 <- lda2\(scaling["LMCD", 1] a2 <- lda2\)scaling[“C12CD”, 1] b <- predict(lda2, newdata = data.frame(LMCD = 0, C12CD = 0))$x - cutoff
Ecuación de la recta: a1LMCD + a2C12CD + b = 0
=> C12CD = -(a1/a2)*LMCD - b/a2
LD1 <- pred2$x[, 1]
ggplot(data.frame(LD1, sexoN = Hombro_sinNA$sexoN),
aes(x = LD1, fill = sexoN)) +
geom_density(alpha = 0.4) +
geom_vline(xintercept = cutoff, linetype = "dashed") +
scale_x_continuous(limits = c(-5, 5)) + # <-- escala de -5 a 5
theme_minimal() +
labs(title = "Distribución sobre la función discriminante",
x = "LD1 (puntuación discriminante)",
y = "Densidad")
Ejemplo de aplicación
Coeficientes del modelo ya calculado
coeficiente LMCD a1 <- -0.0727
coeficiente C12CD a2 <- -0.256
intercepto c <- 18.9578
punto de corte cutoff <- 0.236
Valores nuevos
LMCD_nuevo <- 160
C12CD_nuevo <- 66
Calcular la función discriminante
D_x <- a1 * LMCD_nuevo + a2 * C12CD_nuevo + c
Mostrar resultado
cat("Valor de D(x) para LMCD=160 y C12CD=66:", round(D_x, 4), "\n")
Valor de D(x) para LMCD=160 y C12CD=66: -0.3059
Clasificación
if (D_x > cutoff) {
cat("Clasificación: Hombre\n")
} else {
cat("Clasificación: Mujer\n")
}
4. Modelo logistico para calcular % de clasificación correcta
Omitir datos perdidos
Hombro_sinNA <- na.omit(Hombro[, c("sexoN", "LMCD","C12CD")])
# Ajustar modelo logístico binario
modelo1 <- glm(sexoN ~ LMCD,
data = Hombro_sinNA,
family = binomial(link = "logit"))
summary(modelo1)
# Probabilidades de ser “Hombre” Hombro_sinNA$prob_Hombre2 <- predict(modelo1, type = “response”)
Clasificación usando 0.5 como punto de corte
Hombro_sinNA$predicho2 <- ifelse(Hombro_sinNA$prob_Hombre2 >= 0.5, "Hombre", "Mujer")
# Matriz de confusión
table(Real = Hombro_sinNA$sexoN, Predicho2 = Hombro_sinNA$predicho2)Porcentaje de acierto
mean(Hombro_sinNA$sexoN == Hombro_sinNA$predicho2) * 100
ggplot(Hombro_sinNA, aes(x = prob_Hombre2, fill = sexoN)) + geom_density(alpha = 0.45, linewidth = 1) + geom_vline(xintercept = 0.5, color = “red”, linetype = “dashed”, linewidth = 1) + scale_fill_manual(values = c(“Hombre” = “#1f77b4”, “Mujer” = “#ff7f0e”)) + theme_minimal(base_size = 14) + labs( title = “Distribución de las probabilidades predichas de ser Hombre”, subtitle = “Modelo logístico con LMCD y C12CD”, x = “Probabilidad predicha: P(Hombre)”, y = “Densidad”, fill = “Sexo real” ) + theme( plot.title = element_text(face = “bold”), legend.position = “top” ) coef(modelo1)
Ajustar el modelo logístico
glm_lmcd <- glm(sexoN ~ LMCD,
data = Hombro_sinNA,
family = binomial(link = "logit"))
# Coeficientes
coef_glm_lmcd <- coef(modelo1)
Ejemplo de aplicación (LMCD hipotético, e.g., 194)
LMCD_hipotetico_glm <- 194
logit_ejemplo_glm_lmcd <- coef_glm_lmcd[1] + (coef_glm_lmcd[2] * LMCD_hipotetico_glm)
prob_ejemplo_glm_lmcd <- exp(logit_ejemplo_glm_lmcd) / (1 + exp(logit_ejemplo_glm_lmcd))
clasificacion_glm_lmcd <- ifelse(prob_ejemplo_glm_lmcd >= 0.5, "Hombre", "Mujer")
cat(paste0("**Ejemplo de Aplicación** (LMCD=", LMCD_hipotetico_glm, " mm):\n"))ggplot(Hombro_sinNA, aes(x = prob_Hombre2, fill = sexoN)) + geom_density(alpha = 0.4) + geom_vline(xintercept = 0.5, linetype = “dashed”, color = “black”) + labs(title = “Probabilidades Predichas de ser Hombre (LMCD)”, x = “P(Hombre)”, y = “Densidad”) + theme_minimal()
5. Ajustar modelo logistico, calcular el % de clasificación correcta usando clasificación % correcta usando variables LMCD Y C12CD
Ajustar modelo logístico binario
modelo2 <- glm(sexoN ~ LMCD + C12CD,
data = Hombro_sinNA,
family = binomial(link = "logit"))
summary(modelo2)
Probabilidades de ser “Hombre”
Hombro_sinNA$prob_Hombre <- predict(modelo2, type = "response")
Clasificación usando 0.5 como punto de corte
Hombro_sinNA$predicho <- ifelse(Hombro_sinNA$prob_Hombre >= 0.5, "Hombre", "Mujer")
Matriz de confusión
table(Real = Hombro_sinNA$sexoN, Predicho = Hombro_sinNA$predicho)
Porcentaje de acierto
mean(Hombro_sinNA$sexoN == Hombro_sinNA$predicho) * 100
Gráfica de probabilidades predichas
ggplot(Hombro_sinNA, aes(x = prob_Hombre, fill = sexoN)) +
geom_density(alpha = 0.4) +
geom_vline(xintercept = 0.5, linetype = "dashed") +
theme_minimal() +
labs(title = "Probabilidades predichas de ser Hombre",
x = "P(Hombre)", y = "Densidad")
coef(modelo2)
cat("Ecuación logística:\nlogit(p) = ",
round(coef(modelo2)[1], 4), " + ",
round(coef(modelo2)[2], 4), "*LMCD + ",
round(coef(modelo2)[3], 4), "*C12CD\n")
Ejemplo de aplicación
1. Extraer coeficientes del modelo ya ajustado
b0 <- coef(modelo2)[1]
b1 <- coef(modelo2)[2]
b2 <- coef(modelo2)[3]
2. Valores específicos
LMCD_val <- 160
C12CD_val <- 30
3. Calcular logit y probabilidad
logit_p <- b0 + b1 * LMCD_val + b2 * C12CD_val
probabilidad <- exp(logit_p) / (1 + exp(logit_p))
cat("Probabilidad predicha de ser Hombre =", round(probabilidad, 4), "\n")4. Clasificación según punto de corte = 0.5
clasificacion <- ifelse(probabilidad >= 0.5, "Hombre", "Mujer")
cat("Clasificación predicha =", clasificacion, "\n\n")Clasificación predicha = Mujer
5. Porcentaje de clasificación correcta del modelo completo
Hombro_sinNA$prob_Hombre <- predict(modelo2, type = "response")
Hombro_sinNA$predicho <- ifelse(Hombro_sinNA$prob_Hombre >= 0.5, "Hombre", "Mujer")
porcentaje_clasificacion <- mean(Hombro_sinNA$sexoN == Hombro_sinNA$predicho) * 100
cat("Porcentaje de clasificación correcta del modelo =",
round(porcentaje_clasificacion, 2), "%\n")