1 INTRODUCCIÓN

Una vez organizada la información a través de tablas, el siguiente paso es calcular cierto número de valores característicos que resumen nuestros datos y ayudan a su descripción, ya que también anlizan su distribución en conjunto. Estos valores característicos se calculan a partir de la aplicación de una fórmula matemática.

Estos valores suelen denominarse indistintamente parámetros o estadísticos, pero en realidad no son equivalentes:

Si un estadístico se usa para aproximar un parámetro también se le suele llamar estimador.

Dado que normalmente no es posibles trabajar con toda la población, la expresión habitual asume que trabajamos con estadísticos. Existen diferentes tipos:

Iniciaremos esta actividad estableciendo una semilla para asegurar su reproducibilidad.

set.seed(123)

2 ESTADÍSTICOS DE TENDENCIA CENTRAL

Indican en torno a qué valor (centro) se distribuyen los datos. Son los siguientes:

2.1 LA MEDIA

Cuando hablamos de la media o promedio, normalmente hablamos de la media aritmética, aunque existen otras, como la media geométrica, la armónica etc..

2.1.1 Media aritmética.

Se obtiene sumando todos los datos y dividiéndolos entre el número total de datos:

\[ \overline{x}=\frac{\sum_{i=1}^n x_i}{n}\]

El valor medio de un vector o variable obtiene mediante la función mean().

x <- c(0:10, 50)

mean(x)                        # Con función
sum(x)/length(x)               # Manualmente

La media aritmética tiene las siguientes propiedades:

  • Sólo puede calcularse para variables cuantitativas.

  • Se expresa en las mismas unidades que la variable original.

  • Aunque en su cálculo intervienen todos los valores de la variable original, es sensible a los valores extremos. Para evitar la distorsión estas puntuaciones extremas puedan causar en la media aritmética se utiliza la media recortada (mean α-trimmed), que excluye el α por ciento de los datos. Para ello se añade el argumento trim dentro de la función mean(); en el ejemplo a continuación se excluye el 10 % de los datos del extremo superior y el 10% del extremo inferior.

mean(x, trim = 0.10)

mean(x)
rm(x)

2.1.2 Media aritmética ponderada

Se utiliza cuando no todos los elementos de los que se pretende obtener la media tienen la misma importancia. Su uso es frecuente para promediar porcentajes, tasas, números índices. etc., es decir, cuando la variable presenta variaciones acumulativas. Se calcula multiplicando cada uno de los valores por un peso, que es un valor que indica su importancia con respecto al resto de valores.

Fórmula de la media aritmética ponderada

Ventajas e inconvenientes de la media aritmética ponderada:

  • Cálculo complicado.

  • Los valores extremos tienen menor influencia que en la media aritmética.

  • Cuando en la variable existe un caso con un valor de 0 se anula.

EJEMPLO Para calcular la nota final de una asignatura, se asigna al examen final un valor tres veces superior al resto de las notas. ¿Cuál será la la nota media ponderada de un estudiante que sacó 4,5 y 6,8 en las dos primeras pruebas y 7,4 en el examen final?

notas <- c(4.5,6.8,7.4)
ponderaciones <- c(1,1,3)

El cálculo de la media ponderada se puede realizar manualmente, o bien mediante la función weighted mean().

sum(notas*ponderaciones)/sum(ponderaciones)           # Cálculo manual

weighted.mean(notas, ponderaciones)                   # Función `weighted.mean()`

mean(notas)                                           # Media aritmética 
rm(notas, ponderaciones)

2.1.3 Media geométrica

La media geométrica es más apropiada que la aritmética para calcular el cambio promedio durante el periodo de análisis, típicos del estudio de fenómenos que, como la inflacción o la evolución de precios, están sometido a incrementos/decrecimientos periódicos.

Fórmula de la media aritmética ponderada

EJEMPLO Durante 5 años, el tamaño de un ejample de abejas ha variado con tasas anuales de 14% (primer año), 26% (segundo año etc.), 16%, -38%, -6%. El número de abejas al comienzo del periodo de análisis era 5000. ¿Cuál ha sido el cambio promedio experimentado por el enjambre?.

Si calculamos la media aritmética los resultados serían los siguientes:

abejas <- c(14, 26, 16, -38, -6)              # Tasas de cambio anual en %
mean(abejas)                                  # Media aritmética

De acuerdo con el resultado de la media aritmética, el enjambre habría crecido durante esos 5 años. Para comprobar esta afirmación, procedemos al cálculo de la media geométrica en dos pasos:

  • Paso 1: transformar los porcentajes en tasas de cambio relativo.
crec.relativo <- 1 + abejas/100             
crec.relativo
  • Paso 2: calcular el número de abejas en cada momento multiplicando secuencialmente las tasas de cambio por el número de abejas calculadas para cada momento en el tiempo (eran 5000 al principio).
round(5000 * 
      crec.relativo[1] * 
      crec.relativo[2] * 
      crec.relativo[3] * 
      crec.relativo[4] * 
      crec.relativo[5])

A partir de esos resultados obtenemos una reducción del número de abejas.

También podemos realizar un cálculo análogo a través de la media geométrica

longitud <- length(abejas)                  # Número de observaciones
media.geometrica <- (crec.relativo[1] * 
                       crec.relativo[2] * 
                       crec.relativo[3] * 
                       crec.relativo[4] *
                       crec.relativo[5]) ^ (1/longitud)

round(media.geometrica,3)                   # Resultados redondeados en %

El valor 0.994 es el porcentaje final de abejas respecto al número inicial, que equivale a 4855.

En el caso de R, existe una función geometric.mean() implantada en el paquete psych que proporciona ese mismo porcentaje final.

if (!require("psych")) install.packages("psych")
library("psych")

geometric.mean(crec.relativo)

2.1.4 Media geométrica ponderada

Es la raíz n-ésima del producto de todos los números de un vector o variable. Es recomendada para promediar razones, interés compuesto y números índice.

No existe en R ninguna función que calcule la media geométrica ponderada, pero existe una función ex profeso para su cálculo.Para asegurarnos que la implementación es correcta, retomamos el ejemplo del enjambre de abejas de la sección anterior.

media.geometrica.ponderada <- function(vals, weights){prod(vals)^(1/sum(weights))}

En el ejemplo anterior se observaron durante 5 años las tasas de cambio de la población de abejas. Las tasas de variación anual fueron 1,14, 1,26, 1,16, 0,62, 0,94, que corresponden a x1, x2, … xn. Para reproducir el resultado anterior, establecemos todos los pesos w1, w2, … wn en 1.

round(media.geometrica.ponderada(crec.relativo, rep(1,5)),3)
rm(list = ls())

2.1.5 Media armónica

Es un estadístico propiado para calcular promedios de la velocidad (km/h) o la densidad de población (pob(km2).

Fórmula de la media armónica

EJEMPLO La distancia desde Torrelavega a Santander es 32 km. La velocidad media del vieaje de ida esta mañana ha sido de 105 km/h. A la vuelta, debido a la alta densidad de tráfico, la velocidad era de 60 km/h. ¿Cuál ha sido la velocidad media de todo el viaje?.

ida <- 105
vuelta <- 60

# Media aritmética
mean(c(ida,vuelta))

El resultado de la media aritmética ignora el hecho que se condujo a 105 km/h durante menos tiempo que a 60 km/h (viaje de vuelta). Para calcular la media armónica podemos implementar una función.

media.armonica <- function(vals){length(vals) / sum(1/vals)}
media.armonica(c(ida,vuelta))

O utilizar la función harmonic.mean del paquete psych.

library("psych")
harmonic.mean(c(ida, vuelta))

2.1.6 Media armónica ponderada

Un ejemplo clásico del uso de la media armónica comparada es el cálculo de la densidad de la población (por ejemplo, de un país) según unidades administrativas más pequeñas (por ejemplo, estados, provincias, municipios), con diferente tamaño (población y sobre todo, superficie). El siguiente ejemplo ilustra el problema.

EJEMPLO ¿Cómo calcular la densidad media de la población urbana de las capitales de cada uno de los estados federados de Alemania?.

Para llevar a cabo este cálculo, primero es necesario cargar el fichero desde una página web.

ciudades <- read.csv('https://userpage.fu-berlin.de/soga/200/2010_data_sets/cities.csv')
ciudades
names(ciudades)[1:4] <- c("ciudad", "estado", "area", "poblacion")
  • PASO 1: Creamos una función ad-hoc
media_armonica_ponderada <- function(vals, weights){
  if (sum(weights) != 1){
    weights = weights/sum(weights)
    print('Las ponderaciones no suman 1. Se normaliza las ponderaciones...')
  } 
  sum(weights)/sum(weights/vals)
}
  • PASO 2: Creación de una nueva variable (columna) en la que se calcula la densidad de habitantes (número de habitantes/area en \[km^2\])
ciudades$densidad <- ciudades$poblacion / ciudades$area
ciudades
  • PASO 3: cálculo del peso (ponderación) de cada ciudad, según su población
ciudades$ponderacion <- ciudades$poblacion / sum(ciudades$poblacion)
ciudades
sum(ciudades$ponderacion)             # Comprobación: todas las ponderaciones suman 1
  • PASO 4: cálculo de la media armónica ponderada y comparación con la media aritmética.
media_armonica_ponderada(ciudades$densidad, ciudades$ponderacion)

La densidad media (ponderada) de las capitales de estado en Alemania es 2363 habitantes/km2. La media aritmética sería

mean(ciudades$densidad)
rm(list=ls())

2.2 LA MODA

Es el valor más frecuente en un conjunto de valores, es decir, el que más veces se repite. Si dos puntuaciones adyacentes tienen la frecuencia máxima, se calcula como promedio de esas dos puntuaciones adyacentes; además, si dos o más puntuaciones tienen la misma frecuencia (máxima) hablaremos de una distribución bimodal o multimodal. Este último caso debe analizarse con precaución, ya que puede ser resultado de la mezcla de poblaciones o muestras.Si todas las puntuaciones de un grupo tienen la misma frecuencia no se podría calcular la moda. Se representa por \(Mo\).

Ventajas e inconvenientes:

  • Poco representativa, salvo cuando destaca mucho.

  • No intervienen todos los valores de la distribución.

  • Única medida de posición central que puede obtenerse en las variables de tipo cualitativo.

El vector sobre el que calcularemos la moda es el siguiente:

x <- c(1, 1, 1, 4, 4, 5, 5, 5, 7, 8, 9, 9, 9, 9)

R no dispone de una función en su paquete base que nos permita calcular la moda. Existen, no obstante, algunos paquetes que sí permiten calcular la moda. En este caso instalaremos el paquete modeest y usaremos la función mlv() con el argumento method = "mfv", que devuelve el valor(es) más frecuente en un vector numérico.

if (!require("modeest")) install.packages("modeest")
library(modeest)
mlv(x, 
    method = "mfv")

Otra posibilidad para el cálculo de la moda sería requerir a R qué valor de una tabla tiene una frecuencia máxima.

as.numeric(names(which(table(x)==max(table(x)))))

2.3 LA MEDIANA

Es el valor que ocupa la posición central de la serie de números cuando los datos están ordenados de menor a mayor, y por lo tanto separa una distribución en dos partes iguales. Se representa por \(Me\) y se calcula de la siguiente manera:

  • Si \(n\) es par, la media de los dos datos centrales: \[ \frac{x_{(n/2)}+x_{(n/2+1)}}{2}.\]
  • Si \(n\) es impar, el dato central: \(x_{(n+1)/2}\).
median(x)                                                            
rm(list=ls())

Importa el fichero zonas_verdes.Rdata. A partir de él, calcula los siguientes estadísticos de la(s) variable(s) que consideres más apropiadas.

3 ESTADÍSTICOS DE DISPERSIÓN.

Informan sobre cuánto se alejan los valores del centro de la distribución. Cuanto mayores sean las diferencias entre los individuos del grupo respecto de su promedio, mayor será la dispersión, y por tanto, la utilidad de los valores de tendencia central será menor.

Dispersión alrededor de la media

Los principales estadísticos de dispersión son:

Comenzaremos creando un nuevo vector con los siguientes valores.

x <- c(18,22,16,19,23,18,35,16,45,20,20,22,40,18,45)

3.1 RANGO (RECORRIDO)

Es la diferencia entre el dato de mayor valor y el de menor valor. Por lo tanto, en su cálculo sólo se tiene en cuenta los valores extremos por lo que no necesariamente implica una gran dispersión. En R se calcula anidando la función diff() en la función range().

max(x)-min(x)
range(x)

rango <- diff(range(x))

3.2 DESVIACIÓN RESPECTO A LA MEDIA

Es la diferencia entre cada valor y la media aritmética de toda la distribución.

\[ D_i= x_i - \overline{x}\]

desv_respecto_media <- x - mean(x)

3.3 DESVIACIÓN MEDIA

Es la media aritmética de los valores absolutos de las desviaciones respecto a la media.

Fórmula de la media armónica

desv_media <- sum(abs(x - mean(x)))/length(x)

3.4 VARIANZA

La varianza y su raíz cuadrada positiva, la desviación típica, son las medidas de dispersión más importantes.

La varianza es la media aritmética de las diferencias al cuadrado entre los datos \(x_i\) y la media \(\overline{x}\) de la variable. Se representa por \(σ^2\) (población) o \(s^2\) (muestra).

Fórmula de la varianza

La varianza muestral se calcula con la misma fórmula que la varianza salvo que el denominador es \(n-1\) en lugar de \(n\). \[ \tilde{s}^2 =\frac{\sum_{i=1}^n (x_i-\overline{x})^2}{n-1}. \]

La distinción entre la versión muestral y la “verdadera” de la varianza está motivada por la interrelación entre la estadística descriptiva y la inferencial. Por un lado, se debe mediar la variabilidad de un conjunto de datos cuantitativos mediante su varianza “verdadera”; pero, por otro lado, nuestro conjunto de datos será, normalmente, una muestra de una población mucho mayor, de la que querremos estimar información, y en concreto su variabilidad. Con las técnicas de la estadística descriptiva, resumimos y representamos las características de esta muestra concreta; pero este estudio suele ser solo un paso previo al análisis inferencial, cuyo objetivo no es analizar esta muestra en si misma, sino inferir información sobre toda la población a partir de esta muestra.

La varianza de una muestra tiende a dar valores más pequeños que la varianza real de la población. Para muestras grandes, la diferencia no es sustancial: si \(n\) es grande, dividir por \(n\) o por \(n-1\) no supone una gran diferencia, sobre todo si tenemos en cuenta que se trata de estimar la varianza de la población, no de calcularla exactamente. Pero si el tamaño de la muestra es pequeño (menos de 25 individuos), la varianza muestral de una muestra aproxima significativamente mejor a la varianza real de la población que su varianza “verdadera”.

La varianza tiene las siguientes propiedades matemáticas:

  • Siempre es mayor que 0 \(s^2\geq 0\), porque es una suma de cuadrados de numeros reales.

  • Sólo la varianza equivale a 0 \(s^2=0\) cuando todos los sumandos \((x_i-\overline{x})^2\) son 0 y, por lo tanto, todos los datos son iguales a su media. - Si a todos los individuos de la variable se les suma un número la varianza no varía; si se multiplican por un número la varianza queda multiplicada por el cuadrado de dicho número

  • Es muy sensible a valores atípicos

La función de R para el cálculo de la varianza (muestral) es var().

var(x)

3.5 LA DESVIACIÓN TÍPICA (tipo o estándar)

La desviación típica es la raíz cuadrada \(s\) de la varianza: \(s=\sqrt{s^2}\). Se representa por σ (población) o s (muestra)

La desviación típica muestral es la raíz cuadrada positiva \(\tilde{s}\) de la varianza muestral: \(\tilde{s}=\sqrt{\tilde{s}^2}\).

Fórmula de la desviación típica

Se calcula con la función sd()

sd(x)

Propiedades de la desviación típica:

  • Siempre produce 0 o un valor positivo.

  • Si todos los valores de la variable se les suma un número la desviación típica no varía; si se multiplican por ese número la desviación típica queda multiplicada por dicho número.

Observaciones:

  • Cuanta más pequeña es la desviación típica más concentrados están los datos alrededor de la media.

  • Muy sensible a las puntuaciones extremas.

Calcula los siguientes estadísticos de la variable superficie perteneciente al dataframe zonas_verdes.Rdata. Crea luego una tabla con todos esos valores

  • Rango

  • Desviación media

  • Varianza

  • Desviación típica

3.6 PUNTUACIONES TÍPICAS (estandarizadas)

Se utilizan para comparar variables diferentes, que pueden tener medias y desviaciones típicas muy diferentes.

Fórmula para el cálculo de las puntuaciones típicas Observaciones:

  • La media aritmética de las puntuaciones típicas es 0.

  • La desviación típica de las puntuaciones típicas es 1.

  • Es un estadístico adimensional, o sea, independiente de las unidades utilizadas.

z <- (x - mean(x)) / sd(x)                     
z

EJEMPLO Un problema típico al analizar la evolución de las precipitaciones en España son los valores tan dispares existentes entre unos lugares y otros. Esto supone una dificultad, que se puede solventar transformando las precipitaciones en puntuaciones estandarizadas. Para ilustrar esta solución, analizaremos las evolución de las precipitaciones anuales en 3 observatorios españoles, Santiago de Compostela-Labacolla, Bilbao-Sondica y Murcia-Alcantarilla

url <- "http://personales.unican.es/rasillad/docencia/G14/TEMA_2/precipitaciones.csv"
precipitaciones <- read.csv2(url)

# Cambio en el nombre de la variable
names(precipitaciones)[1] <- "AÑOS"

Los valores medios anuales de precipitación son muy diferentes en estos observatorios

mean(precipitaciones$SANTIAGO)/10
mean(precipitaciones$BILBAO)/10
mean(precipitaciones$MURCIA)/10

Representación gráfica de la evolución de las precipitaciones

plot(precipitaciones$AÑOS, precipitaciones$SANTIAGO/10, lwd = 2, col = "blue", type = "l", lty = 1, ylim = c(0,2700))
lines(precipitaciones$AÑOS, precipitaciones$BILBAO/10, lwd = 2, col = "red", type = "l", lty = 2, add = TRUE)
lines(precipitaciones$AÑOS, precipitaciones$MURCIA/10, lwd = 2, col = "darkgreen", type = "l", lty = 4, add = TRUE)
precipitaciones$ZSANTIAGO <- (precipitaciones$SANTIAGO - mean(precipitaciones$SANTIAGO)) / sd(precipitaciones$SANTIAGO)
precipitaciones$ZBILBAO <- (precipitaciones$BILBAO - mean(precipitaciones$BILBAO)) / sd(precipitaciones$BILBAO)
precipitaciones$ZMURCIA <- (precipitaciones$MURCIA - mean(precipitaciones$MURCIA)) / sd(precipitaciones$MURCIA)

Representación gráfica de la evolución de las precipitaciones transformadas en puntuaciones Z

plot(precipitaciones$AÑOS, precipitaciones$ZSANTIAGO, lwd = 2, col = "blue", type = "l", lty = 1, ylim = c(-3,3))
lines(precipitaciones$AÑOS, precipitaciones$ZBILBAO, lwd = 2, col = "red", type = "l", lty = 2, add = TRUE)
lines(precipitaciones$AÑOS, precipitaciones$ZMURCIA, lwd = 2, col = "darkgreen", type = "l", lty = 4, add = TRUE)

3.7 EL COEFICIENTE DE VARIACIÓN.

Se calcula como la relación entre la desviación típica y la media y se expresa como porcentaje. Permite comparar las dispersiones de dos distribuciones distintas, siempre que sus medias sean positivas. En líneas generales, a mayor coeficiente de variación mayor dispersión

Fórmula del coeficiente de variación

cv <- sd(x) / mean(x) * 100                              
cv

El coeficiente de variación también puede aplicarse a los datos de precipitaciones para saber cuál de los 3 observatorios presenta una mayor regularidad/irregularidad

CV_SANTIAGO <- sd(precipitaciones$SANTIAGO) / mean(precipitaciones$SANTIAGO) * 100                              
CV_BILBAO <- sd(precipitaciones$BILBAO) / mean(precipitaciones$BILBAO) * 100                              
CV_MURCIA <- sd(precipitaciones$MURCIA) / mean(precipitaciones$MURCIA) * 100

De acuerdo con estos datos, el observatorio de Murcia es el que presenta una mayor irregularidad interanual, lo cual es típico de los climas mediterráneos, mientras el de Bilbao es el que presenta una mayor regularidad interanual, típico a su vez de los climas oceánicos.

4 DE POSICIÓN: CUANTILES Y PERCENTILES

Se llama cuantil \(p\) (o de orden \(p\)) (0<\(p\)<100) a aquel valor que, una vez ordenado un conjunto de datos de manera creciente, divide a la variable en \(n\) subconjuntos, dejando por debajo de él un \(p\) por ciento de los datos y por encima el 100−\(p\) por ciento. Por ejemplo, si \(p\)=50, el percentil de orden 50 corresponde a la mediana.

Generalmente, los percentiles van de 1 a 100, mientras los cuantiles se toman de 0 a 1, y es entonces lo mismo el percentil 12, por ejemplo, que el cuantil 0.12.

En R podemos recurrir a su cálculo manual o utilizando varias funciones.

x <- runif(1000,                  # Número de observaciones a ser generadas
           min = 17,              # Límite inferior de la distribución
           max = 89)             # Límite superior de la distribución

EJEMPLO: Cálculo del percentil 38

k <- 25                     # Definimos el percentil 
N <- length(x)              # Número total de datos de la variable
round(k*N/100)              # Posción que corresponde al percentil 38.
sort(x)[round(k*N/100)]     # Valor original que corresponde a esa posición
length(x)
nub.longitud.ordenada <- sort(x)              # Ordenamos los datos

orden25 <- round(25*length(nub.longitud.ordenada)/100)                               
nub.longitud.ordenada[orden25]        # Q25: elemento de una lista ordenada igual o mayor igual que el 25% de los datos

orden50 <- round(50*length(nub.longitud.ordenada)/100)   
nub.longitud.ordenada[orden50]

orden75 <- round(75*length(nub.longitud.ordenada)/100)   
nub.longitud.ordenada[orden75]
quantile(x, probs = 0.38)                   # Cálculo del percentil 38
quantile(x, probs = c(0.25,0.50,0.75))      # Cálculo de los percentiles 25, 50 y 75

También es posible obtener el valor que corresponde a un determinado percentil de un conjunto de datos. Este rango proporciona el porcentaje de valores de ese conjunto de datos cuyo valor es menor que ese citado rango.

Creación de una función

percentile.ranked <- function(a.vector, value) {
  numerator <- length(sort(a.vector)[a.vector < value]) 
  denominator <- length(a.vector)
  round(numerator/denominator,3)*100 
}

Por ejemplo, calculamos el percentil de 63.5

value <- 63.5
percentile.ranked(x, value)

4.1 RANGO O INTERVALO INTERCUARTÍLICO

El rango intercuartílico que es la diferencia entre el cuartil 3º y el cuartil 1º. \(Q_{0.75}-Q_{0.25}\). También se llama a veces rango intercuartílico al intervalo intercuartílico: el intervalo \([Q_{0.25},Q_{0.75}]\).

El IQR es un parámetro muy utilizado para la identificación de datos atípicos. Se define como valor atípico leve aquel que dista 1,5 veces de el rango intercuantílico, bien por debajo de Q1 o bien por encima de Q3

q < Q1 – 1,5 · IQR o bien q > Q3 + 1,5 · IQR

y valor atípico extremo aquel que dista 3 veces el rango intercantílico por debajo de Q1 o por encima de Q3

q < Q1 – 3 · IQR o bien q > Q2 + 3 · IQR

El rango intercuartílico se calcula con la función IQR()

IQR(x)

EJEMPLO Otro problema típico derivado de los valores tan diferentes de precipitación es la definición de sequía y de años anómalos. Ambos fenómenos pueden solventarse mediante el cálculo de los quintiles y del IQR. De nuevo recurrimos a las precipitaciones anuales de Santiago de Compostela-Labacolla, Bilbao-Sondica y Murcia-Alcantarilla.

quantile(precipitaciones$SANTIAGO, probs = c(0.20,0.40,0.60,0.80))
min(precipitaciones$SANTIAGO)
max(precipitaciones$SANTIAGO)
clasificacion_años_SANTIAGO <- cut(precipitaciones$SANTIAGO, 
                                   breaks = c(min(precipitaciones$SANTIAGO),
                                              quantile(precipitaciones$SANTIAGO, probs = c(0.20,0.40,0.60,0.80)),
                                              Inf),  
                                   right = FALSE, 
                                   labels = c("Muy seco", "Seco", "Normal", "Lluvioso", "Muy lluvioso"))
clasificacion_años_SANTIAGO

santiago <- cbind(precipitaciones$SANTIAGO, precipitaciones$ZSANTIAGO, clasificacion_años_SANTIAGO)
santiago

Para identificar si hay años extremos, se puede crear una nueva variable, denominada atipicos recodificando aquellos valores que se cumplan las condiciones señaladas

Q1_SANTIAGO <- quantile(precipitaciones$SANTIAGO, probs = 0.25)
Q3_SANTIAGO <- quantile(precipitaciones$SANTIAGO, probs = 0.75)
IQR_SANTIAGO <- IQR(precipitaciones$SANTIAGO)
atipico_superior_SANTIAGO <- Q3_SANTIAGO + IQR_SANTIAGO * 1.5
atipico_inferior_SANTIAGO <- Q1_SANTIAGO - 1.5 * IQR_SANTIAGO

Finalmente, sólo quedar recodificar los diferentes años

precipitaciones$atipicos_SANTIAGO[precipitaciones$SANTIAGO <= atipico_inferior_SANTIAGO | precipitaciones$SANTIAGO >= atipico_superior_SANTIAGO] <- 1
precipitaciones$atipicos_SANTIAGO[precipitaciones$SANTIAGO > atipico_inferior_SANTIAGO & precipitaciones$SANTIAGO < atipico_superior_SANTIAGO] <- 0
rm(list=ls())

4.1.1 5 números de Tukey

Consiste en cinco números que resumen los estadísticos más importantes para resumir un conjunto de datos, informando sobre su amplitud, su dispersión y el promedio. Además, el resumen de estos cinco números también se puede representar de forma gráfica, lo cual facilita la visualización de estas características de un conjunto de datos,

Los cinco números son:

  • La media.

  • Los cuartiles 1(\(Q_1\)) y 3 (\(Q_3\))

  • Los valores mínimo y máximodel grupo de datos.

fivenum(x)

La función summary() proporciona además la mediana.

summary(x)
rm(list=ls())

5 ESTADÍSTICOS DE FORMA

Muestran si la forma de una distribución presenta ciertas características que permiten clasificarlas en diferentes tipos. Son últiles para compararlas con la forma de distribuciones de probabilidad conocidas, para identificar la que mejor representa el comportamiento de los datos.

Hay dos tipos de estadísticos.

5.1 Estadísticos de asimetría

Son indicadores que permiten establecer el grado de simetría (o asimetría) que presentan los datos de una distribución sin tener que representarlos gráficamente. Para medir la simetría se toma como eje de asimetría la recta paralela al eje de ordenadas que pasa por la media aritmética.

  • Si una distribución es simétrica, existe el mismo número de valores a la derecha que a la izquierda de la media aritmética, por tanto, el mismo número de desviaciones con signo positivo que con signo negativo.

  • La asimetría es positiva (o a la derecha) cuando la “cola” a la derecha de la media es más larga que la de la izquierda, es decir, si hay valores más separados de la media a la derecha.

  • La asimetría es negativa (o a la izquierda) si la “cola” a la izquierda de la media es más larga que la de la derecha, es decir, si hay valores más separados de la media a la izquierda.

Existen diferentes estadísticos para medir la asimetría. Los más importantes son:

  • Coeficiente de asimetría de Fisher.

  • Coeficiente de asimetría de Pearson.

  • Coeficiente de asimetría de Bowley-Yule.

5.1.1 Coeficiente de asimetría de Fisher

El coeficiente de asimetría de Fisher es el más utilizado. Tiene en cuenta el número casos, el promedio y la desviación típica.

Fórmula del coeficiente de asimetría de Fisher

  • Si CAF<0: la distribución tiene una asimetría negativa y se alarga a valores menores que la media.
  • Si CAF=0: la distribución es simétrica.
  • Si CAF>0: la distribución tiene una asimetría positiva y se alarga a valores mayores que la media.

Coeficiente de asimetría de Fisher

Para obtener el coeficiente de asimetría es necesario instalar la librería “moments”. Esta librería contiene la función skewness que proporciona el coeficiente de asimétrica de Fisher.

if (!require("moments")) install.packages("moments")    
library(moments)

x <- sample(17:89,
           1000,                # Número de observaciones a ser generadas
           replace = TRUE)
hist(x)

skewness(x)

5.1.2 Coeficiente de asimetría de Pearson

Es otro de los coeficientes de gran uso, pero limitado a distribuciones uniformes, unimodales y moderadamente asimétricas. Se basa en comparar la media de la distribución con respecto a su moda.

Fórmula del coeficiente de asimetría de Pearson

  • Si CAP<0: la distribución tiene una asimetría negativa, puesto que la media es menor que la moda.

  • Si CAP=0: la distribución es simétrica.

  • Si CAP>0: la distribución tiene una asimetría positiva, ya que la media es mayor que la moda.

Coeficiente de asimetría de Pearson Su cálculo en R es relativamente sencillo, pero necesitamos el cálculo de la moda a través de la función mlv() con el argumento “mfv”.

library(modeest)
cap_pearson <- (mean(x)- mlv(x, method="mfv"))/mean(x)

5.1.3 Coeficiente de asimetría de Bowley-Yule

Está basado en la posición de los cuartiles y la mediana, y utiliza la siguiente expresión:

Coeficiente de asimetría de PearsonD:/G14_EE_CC_SS_2023/TEMA_2_Estadistica_Descriptiva/graficos/3_imagen_Estadisticos_18_Bowley.png){width=40%}

  • Si CAB<0: la distribución tiene una asimetría negativa, puesto que la distancia de la mediana al primer cuartil es mayor que al tercero.

  • Si CAB=0: la distribución es simétrica, ya que el primer y tercer cuartil están a la misma distancia de la mediana.

  • Si CAB>0: la distribución tiene una asimetría positiva, ya que la distancia de la mediana al tercer cuartil es mayor que al primero.

Coeficiente de asimetría de de Bowley

Para obtener el Coeficiente de asimetría de Bowley, es necesario obtener antes los cuartiles.

Q1 <- quantile(x, probs = 0.25)
Q2 <- quantile(x, probs = 0.50)
Q3 <- quantile(x, probs = 0.75)

cab <- (Q3+Q1-2*Q2)/(Q3-Q1)
cab

También existe un Coeficiente Absoluto de Asimetría, que requiere también requiere haber calculado previamente los cuartiles, pero en el que el denominador es la desviación típica.

caa <- (Q3+Q1-2*Q2)/sd(x)
caa

5.2 Estadísticos de apuntamiento (curtosis)

Es una medida de forma que mide cuán escarpada o achatada está una distribución. Indica la cantidad de datos que hay cercanos a la media, de manera que a mayor curtosis, más escarpada (o apuntada) será la curva, al tiempo que existe una relativamente menor frecuencia de valores intermedios.

La fórmula de la curtosis es la siguiente:

Formula de la curtosis

Una distribución puede clasificarse como:

– Leptocúrtica: en la cual los valores están muy agrupados alrededor de la media, por lo que la distribución se presenta bastante apuntada y esbelta.

– Mesocúrtica: posee una concentración moderada de valores alrededor de la media.

– Platicúrtica: está distribución tiene una forma más ancha, pues los valores tienden a estar más dispersos.

Formula de la curtosis

Para calcular curtosis en R y RStudio vamos a usar la función kurtosis.

kurtosis(x)

5.3 Ejemplo desarrollado de asimetría y curtosis en R.

Supongamos el siguiente conjunto de datos:

datos <- c(88, 95, 92, 97, 96, 97, 94, 86, 91, 95, 97, 88, 85, 76, 68)

Para visualizar rápidamente la distribución de valores de estos datos se dibuja un histograma

hist (datos)

El paquete moments ofrece la función jarque.test() , que realiza una prueba de bondad de ajuste que determina si los datos de la muestra tienen asimetría y curtosis que coinciden con una distribución normal. Las hipótesis nula y alternativa de esta prueba son las siguientes:

  • Hipótesis nula: el conjunto de datos tiene una asimetría y una curtosis que coincide con una distribución normal.

  • Hipótesis alternativa: el conjunto de datos tiene una asimetría y una curtosis que no coincide con una distribución normal.

Para realizar la prueba:

jarque.test (datos)

El valor p de la prueba resulta ser 0.05756. Dado que este valor es superior al nivel de significación α = .05, no rechazamos la hipótesis nula. Esto significa que no hay evidencia suficiente para decir que este conjunto de datos tiene una asimetría y una curtosis diferente a la distribución normal.

rm(list=ls())

6 ESTADÍSTICOS DE IGUALDAD (concentración)

Son diferentes indicadores del grado de distribución de la variable y cuantifican el grado de igualdad en el reparto de los valores de una variable. Los más importantes son:

6.1 Índice de Gini

Es un índice desarrollado por Corrado Gini para medir la desigualdad en los ingresos, dentro de un país, pero puede utilizarse para medir cualquier forma de distribución desigual. Es un número entre 0 y 1, donde

  • El valor 0 se corresponde con la perfecta igualdad (todos tienen los mismos ingresos).

  • El valor 1 se corresponde con la perfecta desigualdad (una persona tiene todos los ingresos y los demás ninguno).

Índice de Gini 2015

Supongamos la siguiente distribución de rentas en un país cualquiera. Las rentas han sido agrupadas en una serie de categorías, cuyos límites serían \(L_{i1}-L_i\), \(x_i\) la marca de cada clase, \(n_i\) su frecuencia absoluta y \(N_i\) su frecuencia absoluta acumulada.

Ejemplo desarrollado del índice de Gini

Finalmente,

Ejemplo desarrollado del índice de Gini

Desarrollo del índice de Gini con R.

xi <-c(25,75,125,175,225,275,325,375,425,475)                                   # Vector con valores (en este caso, marcas de clase)
ni <-c(23,72,62,48,19,8,14,7,5,2)                                               # Vector con frecuencias
h <-length(ni) 
h
Ni <-cumsum(ni)
Ni
un <-c(xi*ni)
un
M <-cumsum(un)
M
qi <-c(M/sum(un))*100
qi
pi <-c(Ni/sum(ni))*100
pi

Elaboración de la tabla con todas las variables.

tabla_gini <- cbind(xi,ni,Ni,un,M,pi,qi)
tabla_gini

Aplicación de la fórmula de Gini:

Gini <- sum(pi-qi)/(sum(pi)-100)
Gini

6.2 Curva de Lorenz

Su autor fue Max O. Lorenz en 1905. Es una representación gráfica utilizada para mostrar la distribución relativa de una variable, por ejemplo los ingresos según hogares o personas. La curva representa:

  • Porcentaje acumulado del ingreso (\(pi\)).

  • Porcentaje acumulado de personas u hogares (\(qi\)).

La curva parte del origen (0,0) y termina en el punto (100,100).

  • Si el ingreso estuviera distribuido de manera perfectamente equitativa, la curva coincidiría con la línea de 45 grados que pasa por el origen.

  • Si existiera desigualdad perfecta, o sea, si un hogar o persona poseyera todo el ingreso, la curva coincidiría con el eje horizontal hasta el punto (100,0) donde saltaría el punto (100,100).

En general la curva se encuentra en una situación intermedia entre estos dos extremos.

  • Cuanto más cerca de la diagonal, menor concentración/más homogeneidad en la distribución.

  • Cuanto más cerca de los ejes (parte inferior), mayor concentración/menor homogeneidad

En R podemos dibujarlo de manera manual:

plot(pi,qi, 
     type= "l")
lines(pi, qi, 
      col = 2, 
      lwd = 2, 
      type = "p")
legend("topleft", "curva Lorenz", col = 1:2, lty = 1, box.col = 1)

O también podemos dibujarlo a través de la función curva.lorenz del paquete ineq

library(ineq)
curva.lorenz <- Lc(xi,                   # vector con valores
                   ni)                   # vector con frecuencias
plot(curva.lorenz)

Desafío

Un aspecto frecuentemente analizado por la Geografía Agraria es el tamaño de las explotaciones agrarias. Hay regiones en las que predominan las grandes explotaciones (“latifundios”) mientras que en otras la propiedad está mucho más repartida (“minifundio”). Importa el fichero explotaciones_agrarias.csv. Crea con él un dataframe, calcula el índice de Gini para cada provincia y compáralos entre sí.

datos <- read.csv2(...)
attach(datos)                                    # Función que informa a R que se trabajará con el dataframe datos hasta nueva orden

####################################### Estacion 1
h1 <- length(ni1) 
Ni1 <- cumsum(ni1)
un1 <-c(xi*ni1)
M1 <- cumsum(un1)

qi1 <- c(M1/sum(un1))*100
qi1

pi1 <- c(Ni1/sum(ni1))*100
pi1

# Tabla con todas las variables.
tabla_gini1 <- cbind(xi,ni1,Ni1,un1,M1,pi1,qi1)
tabla_gini1

# Aplicación de la fórmula de Gini:

Gini1 <- sum(pi1-qi1)/(sum(pi1)-100)
Gini1


####################################### Estacion 2
h2 <- length(ni2) 
Ni2 <- cumsum(ni2)
un2 <-c(xi*ni2)
M2 <- cumsum(un2)

qi2 <- c(M2/sum(un2))*100
pi2 <- c(Ni2/sum(ni2))*100

# Tabla con todas las variables.
tabla_gini2 <- cbind(xi,ni2,Ni2,un2,M2,pi2,qi2)
tabla_gini2

# Aplicación de la fórmula de Gini:
Gini2 <- sum(pi2-qi2)/(sum(pi2)-100)
Gini2

####################################### Estacion 3
h3 <- length(ni3) 
Ni3 <- cumsum(ni3)
un3 <-c(xi*ni3)
M3 <- cumsum(un3)

qi3 <- c(M3/sum(un3))*100
pi3 <- c(Ni3/sum(ni3))*100

# Tabla con todas las variables.
tabla_gini3 <- cbind(xi,ni3,Ni3,un3,M3,pi3,qi3)
tabla_gini3

# Aplicación de la fórmula de Gini:
Gini3 <- sum(pi3-qi3)/(sum(pi3)-100)
Gini3
  • ¿Cuál es la provincia en la que domina la pequeña explotación? ¿y cuál tiene un índice de Gini más alto?
rm(list=ls())

7 RESÚMENES ESTADÍSTICOS PARA SUBCONJUNTOS DE DATOS

En muchas ocasiones es necesario calcular un determinado estadístico en función de una o varias variables, con el fin de poder establecer comparaciones entre esos subconjuntos. Por ejemplo, imaginemos la siguiente una base de datos. En primer lugar, generaremos la base de datos con esas variables.

dataset <- data.frame(contaminante = round(rnorm(100, sd = 10, mean = 30)),
                     estacion = sample(c("invierno","primavera","verano", "otoño"), size = 100, replace = TRUE),
                     localidad = sample(paste("localidad", 1:4),
                                        size = 100, 
                                        replace = TRUE))
head(dataset)
str(dataset)

Esta base de datos contiene las siguiente variables:

Existen diferentes posibilidades para calcular esos valores desagregados según otras variables.

7.1 Función aggregate()

La función aggregate() calcula resúmenes estadísticos para subconjuntos de datos, que devuelve en forma de dataframe. Su sintaxis es muy sencilla

Función aggregate

Hay dos posibles ejemplos de uso de aggregate() en función de la salida que proporciona

  • En el caso de utilizar el argumento by, para ver los resultados deberemos crear siempre un objeto, que será un dataframe.

  • En el caso de aplicar una fórmula con la sintaxis variable numérica ~ variable categórica, que evita escribir el nombre completo de las variables, debe incluir a cambio el argumento data. Proporciona en pantalla un dataframe, que podemos o no convertir en objeto.

Por ejemplo, para el cálculo de la concentracion de contaminantes según la estación del año.

promedio_estacional <- aggregate(dataset$contaminante,
                                 by = list(dataset$estacion),     # Variable categórica.
                                 mean)                            # Función
promedio_estacional

colnames(promedio_estacional) <- c("estación del año", "Promedio del contaminante") 
promedio_estacional

La segunda opción nos puede servir para calcular el promedio de la concentración de contaminantes según la localidad.

contaminante_localidad <- aggregate(contaminante ~ localidad,
                                    data = dataset, 
                                    FUN = mean) 
contaminante_localidad

Si quisiéramos calcular los deciles extremos la sintaxis sería la siguiente:

deciles_extremos <- aggregate(contaminante ~ localidad, 
                              data = dataset, 
                              FUN = quantile, 
                              probs = c(0.05, 0.95))

También es posible la agregación simultánea. Hay varias posibilidades:

  • Cuando se desee resumir dos o más variables según un único factor, se utiliza la función cbind para concatenarlas:
aggregate(cbind(dataset$contaminante, dataset$estacion), 
          by = list(dataset$localidad), 
          mean,
          na.rm=TRUE)

aggregate(cbind(contaminante, estacion) ~ localidad,        
          FUN = mean, 
          data = dataset,
          na.rm=TRUE)
  • Cálculo de un estadístico (promedio) correspondiente a una variable según los valores que tomen 2 o más factores.
media_2_factores <- aggregate(dataset$contaminante, 
                              by = list(dataset$estacion, dataset$localidad), 
                              FUN = mean,
                              na.rm=TRUE)
  • Varios estadísticos de una única variable según los valores que tomen 2 o más factores.
diferentes_estadisticos <- aggregate(dataset$contaminante, 
                                     by = list(dataset$estacion, dataset$localidad), 
                                     function(x) c(Suma = sum(x), Media = mean(x), SD = sd(x)))
  • También son posibles los cálculos múlitples mediante la sintaxis de la fórmula.
aggregate(contaminante ~ estacion + localidad, 
          data = dataset, 
          FUN = mean, 
          na.rm=TRUE)

7.2 Función tapply()

Esta función también permite aplicar una función a una matriz, array o data frame calculando estadísticos, pero no sólo por filas, sino también por columas.

Función tapply

La clave consiste en configurar el argumento MARGIN, que puede tomar los valores 1 (por columnas), 2 ( por filas) o c(1, 2) por filas y columnas.

Primero, debemos convertir cada variable en un vector y además transformaremos las columna estación en factor, para que los nombres de los niveles se muestren en la salida de la función cuando se ejecute.

contaminante <- dataset$contaminante
localidad<- dataset$localidad

estacion <- factor(dataset$estacion,
                   levels = c(1, 2, 3, 4),
                   labels = c("invierno", "primavera", "verano", "otoño"))

Precio medio por tipo de producto

contaminación_tapply <- tapply(contaminante, localidad, mean)
contaminación_tapply

Si es necesario se puede acceder a cada elemento de la salida especificando el índice deseado entre corchetes.

contaminación_tapply[2]

El formato de salida puede ser algo diferente (una lista) en caso de establecer el argumento simplify como FALSE.

contaminación_tapply <- tapply(contaminante, localidad, mean, simplify = FALSE)
contaminación_tapply

Y en este caso, para acceder a cualquier elemento de la salida hay que utilizar el signo $ y el nombre del elemento.

contaminación_tapply$`localidad 4`

Finalmente, hay también que tener cuidado con la inclusión de valores NA en nuestros datos. Vamos a simular qué pasa cuando incluimos este tipo de datos. Vamos a crear un nuevo dataset e incluiremos algunos datos ausentes

dataset_NA <- as.data.frame(cbind(contaminante, estacion, localidad))
dataset_NA$contaminante <- as.numeric(dataset_NA$contaminante)
dataset_NA$estacion <- as.integer(dataset_NA$estacion)
str(dataset_NA)

dataset_NA[1, 1] <- NA          # Añadiendo valores NA: fila 1, columna 1.
dataset_NA[2, 3] <- NA          # Añadiendo valores NA: fila 2, columna 3.
dataset_NA

Si queremos calcular el valor medio de la contaminación según estaciones astronómicas

tapply(dataset_NA$contaminante, dataset_NA$estacion, mean)

Dentro de la función tapply se pueden especificar argumentos adicionales de la función que estás aplicando después del argumento FUN. En este caso, la función mean permite especificar el argumento na.rm para que no tome en cuenta esos valores NA. Este argumento por defecto es FALSE.

tapply(dataset_NA$contaminante, dataset_NA$estacion, mean, na.rm = TRUE)

7.3 Función by()

Por último, esta función se puede aplicar a un dataframe separando los resultados según un factor.

contaminante_by <- by(dataset$contaminante,
                      dataset$localidad, 
                      FUN=summary,
                      na.rm = TRUE)                                 
contaminante_by

El mayor inconveniente: hay que transformarlos en un dataframe para exportarlos. Para ello podemos recurrir a la siguiente sintaxis:

contaminante_by <- do.call(rbind, contaminante_by)
contaminante_by

Desafío

A partir del dataframe zonas_verdes calcula:

  • Calcula el número total de parques en cada barrio

  • La superficie media, la superficie mediana, los percentiles 25 y 75, la desviación típica y el coeficiente de variación de la superficie según barrios.

  • El número máximo y mínimo de parques en obras/sin obras según barrio.

8 ESTADÍSTICOS SOBRE DATOS AGRUPADOS

Cuando los ordenadores eran mucho menos potentes que los actuales, era costumbre calcular diferentes estadísticos a partir de datos agrupados si el volumen de datos era muy grande. En realidad, es preferible siempre calcular los estadísticos a partir de los datos sin agrupar, puesto que al agruparlos perdemos información.

No obstante, hay ocasiones en las que los datos se obtienen agrupados. En este tipo de situaciones, se pueden calcular los estadísticos y usarlos como aproximaciones a los que corresponderían a los datos “reales”, que en realidad no conocemos.

La media \(\overline{x}\), la varianza \(s^2\), la varianza muestral \(\tilde{s}^2\), la desviación típica \(s\) y la desviación típica muestral \(\tilde{s}\) de un conjunto de datos agrupados se calculan con las mismas fórmulas que para los datos sin agrupar, excepto que sustituimos cada clase por su marca y la contamos con su frecuencia.

Ejemplo del cálculo de la media, la varianza y la desviación típica sobre datos agrupados

ni <- c(1, 8, 10, 9, 8, 4, 2)
xi <- c(15+10*(0:5), 75)

xini <- as.table(xi * ni)
xi2ni <- as.table(xi^2 * ni)

total <- sum(ni)
total

media <- sum(xini)/total                  
media 

varianza <- sum(xi2ni)/total-media^2      
varianza

desv.tipo <- sqrt(varianza)                                         
desv.tipo

Despejaremos los valores anteriores

rm(desv.tipo, media, varianza, total, xini,xi2ni)

Si las aproximaciones a la mayoría de los estadísticos de tendencia central son relativamente fáciles de calcular, la moda, la mediana y los estadísticos de posición son algo más complicados. Se han planteado diferentes métodos.

Uno de los más sencillos es el siguiente, que requiere inicialmente añadir al dataframe original los límites de clase (al menos, cada límite inferior) y calcular y añadir las frecuencias relativas y frelativas acumuladas, y las frecuencias absolutas acumuladas.

En primer lugar, se reproducen los mismos límites de clase de la tabla utilizada como ejemplo.

a <- 10
L1 <- 10
Linf <- L1 + a*(0:6) 
Lsup <- Linf + a

En segundo lugar, calcularemos las frecuencias relativas y relativas acumuladas, y las frecuencias absolutas acumuladas

Ni <- cumsum(ni)                              
fi <- ni/sum(ni) 
Fi <- cumsum(fi)

Combinamos todos los elementos en una misma tabla en forma de data.frame

datos <- data.frame(Linf, Lsup, xi, ni, fi, Ni, Fi)

Por la misma razón que antes, eliminamos los vectores y nos quedaremos sólo con el dataframe.

rm(a, Linf, Lsup, L1, xi, ni, fi, Ni, Fi)

En lo que se refiere a la moda, se trabajará con el clase o intervalo modal, que es la clase con mayor frecuencia. La fórmula para el cálculo de la moda sobre datos agrupados es la siguiente:

\[ Mo=L_{i-1}+a_i\cdot \frac{n_{i+1}}{(n_{i-1} + n_{i+1})} \] donde:

Por lo tanto, es básico determinar cuál es dicha clase modal. Para ello aplicaremos la siguiente sintaxis:

clase.modal <- datos[which(datos$ni==max(datos$ni)),   ]
clase.modal

A continuación, extraeremos de esa clase modal el límite inferior de la clase y el númeor de casos de esa clase modal.

Li <- clase.modal$Linf                   # límite inferior de la clase modal.
ni <- clase.modal$ni                    # ni número de casos de la clase modal.

Por último, extraeremos el número de casos de las clases inmediatamente anterior y posterior a la modal, así como la amplitud. Para ello, primero debemos saber qué filas de nuestros datos corresponden a esas clases.

datos

ni.menos1 <- datos[2,4]                      # número de casos de la clase inmediatamente inferior a la modal
ni.mas1 <- datos[4,4]                      # número de casos de la clase inmediatamente posterior a la modal
ai  <- datos$Lsup[1] - datos$Linf[1]

El resultado final es:

moda.agrupados <- Li + (ni/(ni.mas1+ni.menos1)) * ai

Tabla para el cálculo de la moda sobre datos agrupados

Para la obtención de la mediana sobre datos agrupados es necesario conocer previamente cuál es la clase mediana, es decir, la clase mediana es la clase o intervalo que contiene el valor del medio de todos los datos ordenados de menor a mayor.

La clase mediana se encuentra en el intervalo cuya frecuencia absoluta acumulada es inmediatamente superior al número obtenido con la siguiente fórmula:

\[ \frac{n+1}{2} \]

Donde \(N\) es el número total de datos.

Otra manera de definir la clase mediana es el primer intervalo donde la frecuencia relativa acumulada \(F_i\) sea mayor o igual a \(0.5\).

La fórmula para obtener una aproximación \(Me\) para la mediana a partir de los datos agrupados es la siguiente:

\[ Me=L_i+A_i\cdot \frac{\frac{N}{2} - N_{i-1}}{n_i} \] En la que:

Iniciamos el cálculo con la obtención de la clase mediana

clase.mediana <- datos[which(datos$Fi>=0.5),   ]
clase.mediana

Dado que el intervalo crítico corresponde a la fila 4, a continuación, extraeremos el límite inferior de esa clase y el número de casos de esa clase mediana.

Li.clase.mediana <- clase.mediana[1,1]                   # límite inferior de la clase modal.
ni.clase.mediana <- clase.mediana[1,4]                    # ni número de casos de la clase modal.

Igualmente, extraeremos el número de casos de las clases inmediatamente anterior y posterior a la modal, así como la amplitud. Para ello, primero debemos saber qué filas de nuestros datos corresponden a esas clases.

Ni.clase.mediana.anterior <- datos[3,6]                      # número de casos de la clase inmediatamente inferior a la modal
n <- sum(datos$ni)
ai  <- datos$Lsup[1] - datos$Linf[1]

El resultado final es:

mediana.agrupados <- Li.clase.mediana + ai * ((N/2 - Ni.clase.mediana.anterior)/n)

8.1 Ejemplo: análisis estadístico de una pirámide de edades

Un ejemplo clásico de información recibida en forma de categorías son los datos de población con los que se elaboran las pirámides de edades. Vamos a usar unos datos agrupados para calcular algunos estadísticos de la distribución de la población española por edades. El fichero original contiene dos columnas: una con los grupos de edad y otra con el volumen de la población. Se importará como un data frame con la función read.csv (también se podría con la función read.table con header=TRUE y sep=",").

datos <- read.csv("https://raw.githubusercontent.com/AprendeR-UIB/AprendeR1/master/datos/cens81.csv", 
                  stringsAsFactors=FALSE,            # Son simplemente etiquetas, por lo que las importamos como caracteres, no como factor.
                  encoding="UTF-8")                  # Garantiza que las "ñ" de las etiquetas ("años") se importan de manera correcta .

str(datos)
datos

Siguiendo el esquema visto anteriormente, cada grupo de edades debe tener asignado un valor numérico como marca de clase.Por ejemplo, la clase “De 0 a 4 años” representa el intervalo de edades [0,5), la clase “De 5 a 9 años” representa el intervalo de edades [5,10), y así sucesivamente; esta amplitud de 5 años se rompe a partir de la clase 85 y más. Por eso, para los grupos de edad inferiores a 84 años, tomaremos como marca de clase el punto medio de sus límites, mientras que para el último, de amplitud indeterminada, tomaremos 90 como marca. Estas marcas de clase se añadirán al data frame como una nueva variable.

datos$marcas <- c(2.5+5*(0:16),90)
head(datos)

Una vez realizado esto, podemos calcular diferents estadísticos:

Total <- sum(datos$Población)                                                # Población total
Total

Edad.media <- sum(datos$Población*datos$marcas)/Total                        # Edad media
Edad.media 

Edad.varianza <- sum(datos$Población*datos$marcas^2)/Total-Edad.media^2      # Varianza
Edad.varianza

Edad.desv.tip <- sqrt(Edad.varianza)                                         # Desviación típica
Edad.desv.tip

Desafío. El fichero piramide_edades.csv recoge el número de habitantes de una serie de países del mundo con diferentes estructuras de edades. Importa dicho fichero, elige un país y contesta a las siguientes preguntas

  • ¿cuál es la edad media de la población de ese país?

  • Qué grupo de población era el más numeroso?

  • Recodifica los grupos de edad en 3 categorías (jóvenes, hasta los 14 años, adultos desde 15 a 64, y ancianos desde 64) y calcula los porcentajes que representan cada uno respecto del total.

Si las aproximaciones a la mayoría de los estadísticos de tendencia central son relativamente fáciles de calcular, la mediana y los de posición son algo más complicados. Se han planteado diferentes métodos. Uno de los más sencillos es el siguiente, que requiere inicialmente calcular las frecuencias absolutas acumuladas, relativas y relativas acumuladas, añadiéndolas al dataframe original.

Int.modal <- datos$Edades[which(datos$Población==max(datos$Población))]      # Intervalo modal
Int.modal
datos$FA.acum <- cumsum(datos$Población)
datos$FR <- round(datos$Población/Total, 3)
datos$FR.acum <- round(datos$FA.acum/Total, 3)
datos

La fórmula para obtener una aproximación \(M\) para la mediana de los datos “reales” a partir de los datos agrupados es la siguiente: \[ M=L_{i}+A_i\cdot \frac{\frac{N}{2}- N_{i-1}}{n_i}. \] En la que

  • \(L_i\) es el intervalo crítico para la mediana, el primer intervalo donde la frecuencia relativa acumulada es igual o mayor que 0.5. En este caso, el intervalo crítico es la clase “De 30 a 34 años”, es decir, \([30,35)\).
  • \(N_{i-1}\), la frecuencia absoluta acumulada del intervalo anterior al crítico (si el intervalo crítico es el primero, tomamos \(N_{i-1}=0\)).
  • \(n_i\), la frecuencia absoluta del intervalo crítico.
  • \(A_i\), su amplitud.
  • \(N\), el número total de datos.

En nuestro ejemplo, \(N=37683361\), \(L_i=30\), \(A_i=5\), \(N_{i-1}=18428647\) y \(n_i=2455314\), por lo que \[ M=30+5\cdot\frac{0.5\cdot 37683361-18428647}{2455314}= 30.8411. \]

Este método también permite aproximar el cuantil \(Q_p\) de los datos “reales” a partir de los datos agrupados con la fórmula siguiente: \[ C_k=L_i+A_i\cdot \frac{\frac{k\cdot N}{4}- N_{i-1}}{n_i}. \]

donde ahora el intervalo crítico es ahora el primer intervalo con frecuencia relativa acumulada mayor o igual que \(p\) y el resto de valores se definen relativos a este intervalo crítico.

De este modo, en nuestro ejemplo, el intervalo crítico para \(Q_{0.25}\) es “De 10 a 14 años”, y en este caso \(L_i}=10\), \(A_i=5\), \(N_{i-1}=6383401\) y \(n_i=3302328\), por lo que \[ Q_{0.25}=10+5\cdot\frac{0.25\cdot 37683361-6383401}{3302328}= 14.6. \] En cuanto al tercer cuantil, \(Q_{0.75}\), el intervalo crítico es “De 50 a 54 años”, por lo que \(L_i=50\), \(A_i=5\), \(N_{i-1}=27547001\) y \(n_i=2265091\), y, por consiguiente, \[ Q_{0.75}=50+5\cdot\frac{0.75\cdot 37683361-27547001}{2265091}=51.58. \]