Procedimientos Estatisticos Comunes

Estadisticas resumen

• Medias y otras funciones resumen

En el siguiente ejemplo se calcula la media de la variable X

mean(datos$X)

## [1] 24.36364

La función mean() acepta como argumentos un vector numerico o un dataframe numerico (objetos de fecha tambien son aceptados).

ó invocando a la libreria mosaic, se pueden calcular los estadísticos descriptivos de la variable, es decir, valores mínimos y máximos, cuartiles y número de observaciones

favstats(X, data=datos)

##  min  Q1 median   Q3 max     mean      sd  n missing
##    1 6.5     22 44.5  49 24.36364 19.2472 11       0

Además se podría calcular un resumen de la data completa, usando el comendo summary, como se hace a continuación

summary(datos)

##       caso            X               Y        
##  Min.   : 1.0   Min.   : 1.00   Min.   : 4.00  
##  1st Qu.: 3.5   1st Qu.: 6.50   1st Qu.: 7.00  
##  Median : 6.0   Median :22.00   Median : 8.00  
##  Mean   : 6.0   Mean   :24.36   Mean   :13.36  
##  3rd Qu.: 8.5   3rd Qu.:44.50   3rd Qu.:24.00  
##  Max.   :11.0   Max.   :49.00   Max.   :26.00

-Funciones Similares:

Acá se usa la instrucción attach para no tener que estar llamando a la data cada vez que se necesite calcular algún estadístico a una variable

attach(datos)

median(X)

## [1] 22

var(X)

## [1] 370.4545

sd(X)

## [1] 19.2472

min(X)

## [1] 1

-

max(X)

## [1] 49

sum(X)

## [1] 268

prod(X)

## [1] 2.969091e+12

range(X)

## [1]  1 49

Medias Ponderadas y otras estadisticas

wtd.mean(X,weights = NULL) # calcula la media simple

## [1] 24.36364

wtd.mean(X,weights = c(1,2,1,4,2,6,7,8,9,2,11)) #calcula la media ponderada

## [1] 33.41509

-Verifiquemos este resultado creando nuestro propio código

weights = c(1,2,1,4,2,6,7,8,9,2,11)
suma=0
i=0
for(i in 1:length(X)){
 suma=suma+X[i]*weights[i]
 d=cumsum(weights)
 n=d[length(X)]
 }
    i=i+1
mediaPonderada=suma/n
mediaPonderada

## [1] 33.41509

-

wtd.var(X,weights = c(1,2,1,4,2,6,7,8,9,2,11))

## [1] 255.5552

wtd.quantile(X,weights = c(1,2,1,4,2,6,7,8,9,2,11))

##   0%  25%  50%  75% 100% 
##    1   21   44   45   49

wtd.rank(X,weights = c(1,2,1,4,2,6,7,8,9,2,11))

##  [1]  1.0  2.5  4.0  6.5 15.5 11.5 20.0 27.5 36.0 41.5 48.0

La función wtd.mean() en el paquete Hmisc calcula la media ponderada. Otra de las funciones asociadas son: wtd.var(), wtd.quantile() y wtd.rank()

Momentos

Mientras que los coeficientes de asimetría y curtosis son menos utilizados que la media y la deviación estandar, en ciertos momentos pueden ser útiles.

El coeficiente de asimetría es descrito como el tercer momento alrededor de la media y nos informa si la distribución es simetrica (coeficiente=0).

Curtosis es una función de momento central cuatro.

skewness(X) #coeficiente de asimetria

## [1] 0.1149118

kurtosis(X)

## [1] 1.450516

all.moments(X)

## [1]   1.00000  24.36364 930.36364

all.moments(X, order.max = 5,central = TRUE) 

## [1] 1.000000e+00 3.229740e-16 3.367769e+02 7.101953e+02 1.645156e+05
## [6] 4.997955e+05

El paquete moments facilta el calculo de los coeficientes de asimetría y curtosis dentro de R, así como tambien los momentos de mayor orden, utilizando all.moments()

Media truncada

A veces la data tiene valores muy pequeños o muy grandes con respecto a todos los demás, y para eliminar la atracción que siente la media por valores atípicos, se trunca (elimina) un porcentaje de ellos

mean(X,trim = 0.2) #eliminó 20%

## [1] 24

mean(X,trim=0.5)#Eliminó 50%

## [1] 22

El valor trim puede tomar un rango de 0 a 0.5 y especifica la fracción de observaciones que deben ser cortadas de cada lado de X antes de que la media sea calculada.

frac=0.5 da como resultado la mediana.

Cuantiles

quantile(X, c(.025, .975)) # Calcula los cuantiles del 25% y 97,5% respectivamente

##  2.5% 97.5% 
##  1.25 48.75

quantile(X, seq(from=.95, to=.975, by=.0125)) # Calcula los cuantiles desde el 95% hasta el 97.5% en pasos de 1,25

##    95% 96.25%  97.5% 
## 48.500 48.625 48.750

-Centrando, Normalizando

(X-mean(X)/sd(X)) #normalizando

##  [1] -0.2658278  0.7341722  1.7341722  8.7341722 20.7341722 19.7341722
##  [7] 21.7341722 42.7341722 43.7341722 46.7341722 47.7341722

La función scale() puede operar en matrices y dataframes.

Media e intervalo de confianza 95%

Función cuantil de una distribución t

qt(p, df, ncp, lower.tail = TRUE, log.p = FALSE)

p: vector de probabilidades

df:grados de libertad

lower.tail: cola inferior; si es TRUE (por defecto), las probabilidades son P [X ??? x], de lo contrario, P [X> x].

-

tcrit = qt(.975, df=length(X)-1)
tcrit

## [1] 2.228139

mean(X) + c(-tcrit, tcrit)*sd(X)/sqrt(length(X))

## [1] 11.43319 37.29408

o

t.test(X)$conf.int

## [1] 11.43319 37.29408
## attr(,"conf.level")
## [1] 0.95

Mientras que el intervalo de confianza 95% puede ser generado en términos de la media y la desviación estandar, es mucho mas directo aplicar la función t.test().

Proporciones e intervalo de confianza 95%

heads <- rbinom(1, size = 100, prob = .5)
heads

## [1] 45

binom.test(heads,100)

## 
## 
## 
## data:  heads out of 100
## number of successes = 45, number of trials = 100, p-value = 0.3682
## alternative hypothesis: true probability of success is not equal to 0.5
## 95 percent confidence interval:
##  0.3503202 0.5527198
## sample estimates:
## probability of success 
##                   0.45

prop.test(heads, 100)

## 
##  1-sample proportions test with continuity correction
## 
## data:  heads out of 100
## X-squared = 0.81, df = 1, p-value = 0.3681
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
##  0.3514281 0.5524574
## sample estimates:
##    p 
## 0.45

La función binom.test() calcula un intervalo de confianza Clopper-Pearson exacto basado en la distribución F usando el primer argumento como el numero de exitos y el segundo argumento como el número de eventos.

La función prop.test() calcula un intervalo de confianza aproximado.

La opción conf.level puede ser usada para cambiar el nivel de confianza.

Estimación de parametros

Estadisticas bivariantes

-Correlación

pearsoncorr = cor(X, Y)
pearsoncorr

## [1] -0.1923683

spearmancorr = cor(X, Y, method="spearman")
spearmancorr

## [1] -0.03218603

kendalltau = cor(X, Y, method="kendall")
kendalltau

## [1] 0.01888137

Especificando method=“spearman” o method=“kendall” como una opción para cor() genera los coeficientes de correlación de Spearman o Kendall respectivamente.

Una matriz de variables puede ser utilizada para generar la correlación entre el set de variables.

Tablas de Contingencia

Las Tablas de Contingencia muestra los grupos de miembros a traves de variables categoricas (agrupando)

data(infert, package = "datasets")
CrossTable(infert$education, infert$induced, expected = TRUE)

## Warning in chisq.test(t, correct = FALSE, ...): Chi-squared approximation
## may be incorrect

## 
##  
##    Cell Contents
## |-------------------------|
## |                       N |
## |              Expected N |
## | Chi-square contribution |
## |           N / Row Total |
## |           N / Col Total |
## |         N / Table Total |
## |-------------------------|
## 
##  
## Total Observations in Table:  248 
## 
##  
##                  | infert$induced 
## infert$education |         0 |         1 |         2 | Row Total | 
## -----------------|-----------|-----------|-----------|-----------|
##           0-5yrs |         4 |         2 |         6 |        12 | 
##                  |     6.919 |     3.290 |     1.790 |           | 
##                  |     1.232 |     0.506 |     9.898 |           | 
##                  |     0.333 |     0.167 |     0.500 |     0.048 | 
##                  |     0.028 |     0.029 |     0.162 |           | 
##                  |     0.016 |     0.008 |     0.024 |           | 
## -----------------|-----------|-----------|-----------|-----------|
##          6-11yrs |        78 |        27 |        15 |       120 | 
##                  |    69.194 |    32.903 |    17.903 |           | 
##                  |     1.121 |     1.059 |     0.471 |           | 
##                  |     0.650 |     0.225 |     0.125 |     0.484 | 
##                  |     0.545 |     0.397 |     0.405 |           | 
##                  |     0.315 |     0.109 |     0.060 |           | 
## -----------------|-----------|-----------|-----------|-----------|
##          12+ yrs |        61 |        39 |        16 |       116 | 
##                  |    66.887 |    31.806 |    17.306 |           | 
##                  |     0.518 |     1.627 |     0.099 |           | 
##                  |     0.526 |     0.336 |     0.138 |     0.468 | 
##                  |     0.427 |     0.574 |     0.432 |           | 
##                  |     0.246 |     0.157 |     0.065 |           | 
## -----------------|-----------|-----------|-----------|-----------|
##     Column Total |       143 |        68 |        37 |       248 | 
##                  |     0.577 |     0.274 |     0.149 |           | 
## -----------------|-----------|-----------|-----------|-----------|
## 
##  
## Statistics for All Table Factors
## 
## 
## Pearson's Chi-squared test 
## ------------------------------------------------------------
## Chi^2 =  16.53059     d.f. =  4     p =  0.002383898 
## 
## 
## 

library(mosaic)

-Tabulate categorical data

tally(~ Y + X, margins=TRUE, data=datos)

##        X
## Y        1  2  3 10 21 22 23 44 45 48 49 Total
##   4      0  1  0  0  0  0  0  0  0  0  0     1
##   5      0  0  0  0  0  0  0  0  1  0  0     1
##   6      0  0  0  0  0  0  0  0  0  0  1     1
##   8      1  0  1  0  0  0  0  0  0  1  0     3
##   9      0  0  0  0  0  0  0  1  0  0  0     1
##   23     0  0  0  1  0  0  0  0  0  0  0     1
##   25     0  0  0  0  1  1  0  0  0  0  0     2
##   26     0  0  0  0  0  0  1  0  0  0  0     1
##   Total  1  1  1  1  1  1  1  1  1  1  1    11

La función CrossTable() in el paquete gmodels provee una manera fácil de generar tablas cruzadas.

Mostrando categoria de valores perdidos en una tabla

Puede ser útil mostrar una tabla donde se incluye a los valores pérdidos como una categoria separada.

a <- rep(c(NA, 1/0:3), 10)
a

##  [1]        NA       Inf 1.0000000 0.5000000 0.3333333        NA       Inf
##  [8] 1.0000000 0.5000000 0.3333333        NA       Inf 1.0000000 0.5000000
## [15] 0.3333333        NA       Inf 1.0000000 0.5000000 0.3333333        NA
## [22]       Inf 1.0000000 0.5000000 0.3333333        NA       Inf 1.0000000
## [29] 0.5000000 0.3333333        NA       Inf 1.0000000 0.5000000 0.3333333
## [36]        NA       Inf 1.0000000 0.5000000 0.3333333        NA       Inf
## [43] 1.0000000 0.5000000 0.3333333        NA       Inf 1.0000000 0.5000000
## [50] 0.3333333

table(a,useNA = "ifany")

## a
## 0.333333333333333               0.5                 1               Inf 
##                10                10                10                10 
##              <NA> 
##                10