Practica 2
Nohemi Angel
2025-12-06
INTRODUCCION
En esta practica se realizara comparaciones, construcciones de funcion discriminante, ajustes de modelos logisticos y aplicacion para valores hipoteticos con dos varialbes.
- Abriendo paquetes estableciendo directorio de trabajo y leyendo archivo
setwd("~/Blanca")
##Abriendo paquete pacman
library(pacman)
## El archivo esta en formato SPSS, lo abrimos mediante la libreria haven
p_load(haven,dplyr,ggplot2,MASS)
p_load(tinytex)
Hombro <- read_sav("Datos hombro.sav")## Definimos como factor la variable sexoN
Hombro$sexoN <- factor(Hombro$sexoN,
levels = c(2, 1),
labels = c("Mujer", "Hombre"))
table (Hombro$sexoN) ## Frecuencias de sexo##
## Mujer Hombre
## 30 50Podemos observar que la muestra consta de 80 individuos de los cuales 50 son hombres y 30 mujeres
Estadística descriptiva
Empezemos con la estadistica descriptiva, para resumir por sexo la información de Longitud Máxima de la clavicula derecho
res_lmcd <- Hombro %>%
group_by(sexoN) %>%
summarise(
n = sum(!is.na(LMCD)),
media = mean(LMCD, na.rm = TRUE),
sd = sd(LMHD, na.rm = TRUE)
) %>%
mutate(across(c(media, sd), ~round(.x, 2)))
res_lmcd## # A tibble: 2 × 4
## sexoN n media sd
## <fct> <int> <dbl> <dbl>
## 1 Mujer 25 131. 15.7
## 2 Hombre 35 151. 13.9De la tabla anterior se muestra diferencias significativas entre hombres y mujeres, en las que se observan numeros ligeramente menores que el total debido a datos faltantes, con 25 mujeres y 35 hombres, donde la media en mujeres es aproximadamente de 131 mm, mientras que en hombres es de 151 mm,por lo que podemos analizar que los hombres presentan claviculas mas largas.La desviacion estandar es muy similar en ambos grupos, lo cual indica en la tabla un marcado diformismo sexual en el tamaño de la clavicula.
Ahora veremos el caso de la circunferencia del punto medio de la clavicula derecho
res_cmcd <- Hombro %>%
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_cmcd## # A tibble: 2 × 4
## sexoN n media sd
## <fct> <int> <dbl> <dbl>
## 1 Mujer 29 30.3 2.83
## 2 Hombre 46 35.9 2.33Observamos que en este caso se incluyen 29 mujeres y 46 hombres, en la que la media de las mujeres es de 30.3 mm y en hombres es de 35.6 mm, mostrando que los hombres tienen valores mas grandes.Asimismo podemos observar que la desviacion es mas baja en hombres,lo que podria confirmar que la circunferencia clavicular presenta diformismo sexual evidentemente marcado.
Realizaremos una gráfica comparativa de Longitud máxima del húmero derecho por Sexo
ggplot(Hombro, aes(x = LMCD, fill = sexoN)) +
geom_density(alpha = 0.5) +
labs(
title = "Gráfica 1. Longitud máxima de la clavicula derecho por sexo",
x = "Longitud máxima de la clavicula derecho (mm)",
y = "Densidad"
) +
theme_minimal()## Warning: Removed 20 rows containing non-finite outside the scale range
## (`stat_density()`).
Otra gráfica, ahora de caja.
ggplot(Hombro, aes(x = sexoN, y = LMCD, fill = sexoN)) +
geom_boxplot(alpha = 0.7) +
labs(
title = "Longitud máxima de la clavicula derecho",
x = "Sexo",
y = " "
) +
theme_minimal()## Warning: Removed 20 rows containing non-finite outside the scale range
## (`stat_boxplot()`).
Estas dos graficas nos permiten visualizar de manera clara las diferencias encontradas en los anteriores analisis, observamos por ejemplo en la grafica de densidad, nos indica que la distribucion de los hombres se desplaza hacia valores mas altos que en las mujeres y la grafica de caja muestra que las medidas son mayores en hombres. Esto confirma vusualmente el diformismo sexual de las variables estudiadas. # Análisis discriminante
En este apartado aplicaremos la técnica de análisis discriminante que nos permitirá clasificar un individuo como hombre o mujer dependiendo de los valores discriminantes obtenidos. Un paso intermedio, después de estimar la función discriminante es obtener el punto de corte. Pero para dos variables.
##Ajustar Disciminante Lineal
# Eliminar NAs de forma explícita
Hombro_sinNA <- 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 = Hombro_sinNA)
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[Hombro_sinNA$sexoN == "Hombre"])
centroide_M <- mean(pred$x[Hombro_sinNA$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")## Punto de corte= -0.1791 Si D> -0.1791 es Hombre # Tabla cruzada
tabla_clas <- table(Observado = Hombro_sinNA$sexoN,
Predicho = pred$class)
tabla_clas## Predicho
## Observado Mujer Hombre
## Mujer 22 3
## Hombre 3 32 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:Supongamos que se encuentra una clavicula derecha cuya longitud es de 120.5 mm. usaremos el modelo para predecir
prede1 <- predict(lda1,newdata = data.frame(LMCD = 120.5))
prede1$x## LD1
## 1 -2.430437if(prede1$x>cutoff) cat("Debido a que", prede1$x,">",cutoff, "es hombre") else cat("Debido a que", prede1$x,"<",cutoff,"es mujer")## Debido a que -2.430437 < -0.1790721 es mujerAhora ajustaremos la discriminacion lineal con dos variables
##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=a*LMCD+b*C12CD+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## Call:
## lda(sexoN ~ LMCD + C12CD, data = Hombro_sinNA)
##
## Prior probabilities of groups:
## Mujer Hombre
## 0.4166667 0.5833333
##
## Group means:
## LMCD C12CD
## Mujer 131.3026 30.31206
## Hombre 151.0263 35.77151
##
## Coefficients of linear discriminants:
## LD1
## LMCD 0.0727066
## C12CD 0.2559868 # Predicción de valores mediante modelo
pred2 <- predict(lda2)
# Matriz de confusión
table(Real = Hombro_sinNA$sexoN, Predicho = pred2$class)## Predicho
## Real Mujer Hombre
## Mujer 24 1
## Hombre 1 34 # Porcentaje de acierto
mean(pred2$class == Hombro_sinNA$sexoN)*100## [1] 96.66667 # valor de D cuando LMCD=0 y C12CD=0 (intercepto)
b <- predict(lda2,
newdata = data.frame(LMCD = 0, C12CD = 0))$x
b## LD1
## 1 -18.95782 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## [1] -0.2359658 cat("Función discriminante: D(x) = ", round(a1, 4), " * LMCD + ", round(a2, 4),"*C12CD +",round(b, 4), "\n")## Función discriminante: D(x) = 0.0727 * LMCD + 0.256 *C12CD + -18.9578 cat("Punto de corte=", round(cutoff, 4), "Si D>",round(cutoff, 4), "es Hombre")## 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: LMCD vs C12CD")
a1 <- lda2$scaling["LMCD", 1]
a2 <- lda2$scaling["C12CD", 1]
b <- predict(lda2, newdata = data.frame(LMCD = 0, C12CD = 0))$x - cutoffEjemplo: Supongamos que se encuentra una clavicula derecha cuya longitud maxima es de 122 y una circunferencia del punto medio de la clavicula derecha es de 11. Usaremos el modelo para precedir
ejemplo2 <- predict(lda2, newdata = data.frame(LMCD = 122, C12CD = 11))
ejemplo2$x## LD1
## 1 -7.271756if(ejemplo2$x > cutoff){
cat("El individuo sería clasificado como Hombre")
} else{
cat("El individuo sería clasificado como Mujer")
}## El individuo sería clasificado como MujerRegresión Logística
Aplicaremos la regresión logística para determinar la probabilidad de hombre o mujer usando los huesos del hombro (LMCD Y C12CD). Aplicaremos primero la regresión con una variable (LMCD)
### Ahora aplicaremos Regresión Logística
## Omitir datos perdidos
Hombro_sinNA <- na.omit(Hombro[, c("sexoN", "LMCD")])
# Ajustar modelo logístico binario
modelo1 <- glm(sexoN ~ LMCD,
data = Hombro_sinNA,
family = binomial(link = "logit"))
summary(modelo1)##
## Call:
## glm(formula = sexoN ~ LMCD, family = binomial(link = "logit"),
## data = Hombro_sinNA)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -41.72697 10.79320 -3.866 0.000111 ***
## LMCD 0.29945 0.07717 3.880 0.000104 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 81.503 on 59 degrees of freedom
## Residual deviance: 32.430 on 58 degrees of freedom
## AIC: 36.43
##
## Number of Fisher Scoring iterations: 6 # Probabilidades de ser "Hombre"
Hombro_sinNA$prob_Hombre <- predict(modelo1, 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)## Predicho
## Real Hombre Mujer
## Mujer 3 22
## Hombre 32 3 # Porcentaje de acierto
mean(Hombro_sinNA$sexoN == Hombro_sinNA$predicho) * 100## [1] 90 ## 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(modelo1)## (Intercept) LMCD
## -41.7269684 0.2994501 cat("Ecuación logística:\nlogit(p) = ",
round(coef(modelo1)[1], 4), " + ",
round(coef(modelo1)[2], 4), "*LMCD\n")## Ecuación logística:
## logit(p) = -41.727 + 0.2995 *LMCDAquí vamos a aplicar la regresión lineal con dos variables
### Ahora aplicaremos Regresión Logística
## Omitir datos perdidos
Hombro_sinNA <- na.omit(Hombro[, c("sexoN", "LMCD","C12CD")])
# Ajustar modelo logístico binario
modelo1 <- glm(sexoN ~ LMCD + C12CD,
data = Hombro_sinNA,
family = binomial(link = "logit"))
summary(modelo1)##
## Call:
## glm(formula = sexoN ~ LMCD + C12CD, family = binomial(link = "logit"),
## data = Hombro_sinNA)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -69.0885 22.9591 -3.009 0.00262 **
## LMCD 0.3067 0.1117 2.747 0.00602 **
## C12CD 0.8022 0.3405 2.356 0.01847 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 81.503 on 59 degrees of freedom
## Residual deviance: 18.731 on 57 degrees of freedom
## AIC: 24.731
##
## Number of Fisher Scoring iterations: 8 # Probabilidades de ser "Hombre"
Hombro_sinNA$prob_Hombre <- predict(modelo1, 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)## Predicho
## Real Hombre Mujer
## Mujer 1 24
## Hombre 34 1 # Porcentaje de acierto
mean(Hombro_sinNA$sexoN == Hombro_sinNA$predicho) * 100## [1] 96.66667 ## 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(modelo1)## (Intercept) LMCD C12CD
## -69.0884686 0.3067252 0.8022468 cat("Ecuación logística:\nlogit(p) = ",
round(coef(modelo1)[1], 4), " + ",
round(coef(modelo1)[2], 4), "*LMCD + ",
round(coef(modelo1)[3], 4), "*C12CD\n")## Ecuación logística:
## logit(p) = -69.0885 + 0.3067 *LMCD + 0.8022 *C12CDEjemplo: Supongamos que se encuentra una clavicula derecha cuya longitud maxima es de 120 y una circunferencia del punto medio de la clavicula derecha es de 15. Usaremos el modelo para precedir
nuevo <- data.frame(LMCD = 120, C12CD = 15)
prob_nuevo <- predict(modelo1, newdata = nuevo, type = "response")
cat("Probabilidad de ser Hombre =", round(prob_nuevo, 4), "\n")## Probabilidad de ser Hombre = 0if(prob_nuevo >= 0.5){
cat("El individuo sería clasificado como Hombre")
} else{
cat("El individuo sería clasificado como Mujer")
}## El individuo sería clasificado como MujerEn conclusion, en la presente practica se aplicaron diversos metodos estadisticos para poder analizar el diformismo sexual a partir de dos mediciones osteometricas, se tomaron en cuenta la variable LMCD y C12CD.
El analisis con la variable LMCD alcanzo asta un 90%, mientras que al incluir las dos variables el porcentaje aumento un 96.6%, lo que evidencia que la combinacion de medidas mejoro la presicion del metodo, lo mismo ocurrio con el modelo de regresion logistica mostrando coeficcientes significativos con una presicion del 96.6%
Esto demuestra que las dos variables son utiles y consistentes mediante tecnicas estadisticas.