Modelos Estadísticos 1

Quick-R

Posted in Uncategorized by hector on junio 1, 2009

Quick-R es un website con breve referencias y ejemplos para utilizar R. Denle un vistazo.

Tagged with:

Proyecto final, lineamientos

Posted in Uncategorized by hector on mayo 26, 2009

Estos son los lineamientos que deben seguir en relación al proyecto final.


Objetivos

  • Los estudiantes deberán emplear los conocimientos adquiridos para resolver problemas de análisis de datos.
  • Los estudiantes tendrán que plantear un problema, implementar su solución, y presentarla en forma clara y ordenada en un reporte final.

Requisitos del reporte final

  • Carátula en la que indiquen: nombre de la institución (e.g. UVG), nombre del curso, nombre del estudiante (con carnet y carrera), título del proyecto, y fecha de entrega.
  • Introducción: acá se describe el problema a tratar y algunos de los aspectos importantes relacionados al mismo.
  • Descripción de los datos: acá se especifican claramente las fuentes de la información; qué tipo de datos (numéricos, categóricos) se están manipulando y por qué es el tipo adecuado; qué instrumentos fueron empleados para recopilar los datos (encuestas, observaciones, etc.) ésto debe quedar indicado aunque no sean ustedes quienes recopilen los datos.
  • Estadísticas descriptivas: acá generan tablas y gráficas que ayuden a entender la relación entre las variables estudiadas. Nota: No deben insertar tablas y gráficas que no tengan una descripción textual adecuada. La interpretación de las tablas y gráficas debe estar incluida en el reporte. También deben incluir el código de R con el cual generaron estos resultados. Sean juiciosos, no incluyan tablas y/o gráficas que no son relevantes.
  • Análisis de datos: acá se especifica el tipo de análisis empleado. Deben enunciar claramente las hipótesis (nula y alternativa) e interpretar adecuadamente la(s) prueba(s) realizada(s), i.e. se acepta la hipótesis nula, etc.
  • Conclusión: acá se presenta una sumaria de los resultados. En retrospectiva, deben reflexionar y estipular posibles fallas en el estudio así como proponer mejores.
  • Bibliografía: acá se listan todas las referencias (libros, websites, etc.) que emplearon durante el desarrollo de esta investigación.

En general, traten de seguir los lineamientos de la Guía para la redacción de tesis e informes académicos (UVG, 2001).

Consideren que se tomará muy en cuenta la presentación del reporte.

Fecha de entrega

Jueves, 4 de Junio.


Conjuntos de datos

Pueden emplear datos recopilados por ustedes, o bien, datos disponibles en internet. De cualquier forma debe quedar clara la procedencia de los mismos. Algunas páginas en las que pueden encontrar datasets:

Asegúrense que el dataset que seleccionen contenga al menos 30 observaciones (usualmente filas) y 3 variables (usualmente columnas).

Artículo de R en el New York Times

Posted in Uncategorized by hector on mayo 22, 2009

Les recomiendo este artículo acerca de R que apareció en la sección de tecnología del New York Times:

R is [...] a popular programming language used by a growing number of data analysts inside corporations and academia. [...] Companies as diverse as Google, Pfizer, Merck, Bank of America, the InterContinental Hotels Group and Shell use it.

Tagged with:

Obteniendo resultados estadísticos con R, parte IV

Posted in Tutoriales, Uncategorized by hector on mayo 19, 2009

Esta es la cuarta parte de una breve introducción a los comandos del sistema estadístico R. En esta parte trabajaremos con algunas ideas fundamentales del análisis estadístico y de prueba de hipótesis.


Regresión y Correlación

En muchas ocasiones tendremos variables (numéricas) que están relacionadas de manera lineal. Además, es importante notar que si las variables no tienen una relación lineal, dicha relación puede asumirse en una región de interés (i.e. se toma una aproximación lineal). El modelo de regresión lineal está dado por:
\displaystyle y_i = \alpha + \beta x_i + \varepsilon_i
donde 1 \leq i \leq N. Noten que x_i, y_i representan los distintos valores de las variables x,y; \varepsilon_i representa los errores asociados con dichas mediciones. Estos errores se asumen independientes y obtenidos de una distribución normal N(0,\sigma^2). Asociando esta expresión con la ecuación punto-pendiente de una recta y = m x + b, tenemos que \beta es la pendiente y se denota coeficiente de regresión; además \alpha es el intersecto. Nota: los parámetros \alpha, \beta, \sigma pueden estimarse empleando el Método de Mínimos Cuadrados en donde el objetivo es minimizar el cuadrado de los residuos:
\displaystyle E = \sum_{i=1}^{N} \left(y_i - (\alpha + \beta x_i)\right)^2
después de un poco de Cálculo llegamos a las expresiones:
\hat \beta = \frac{\sum (x_i - \bar x)(y_i - \bar y)}{\sum (x_i - \bar x)^2}
\hat \alpha = \bar y - \hat \beta \bar x

Un concepto bastante relacionado a la idea de regresión lineal es el de correlación lineal. Un coeficiente de correlación es una medida de la asociación entre dos variables; puede tomar valores entre -1 y +1 (donde 0 significa que no hay correlación entre las variables). Nota: “Correlación no implica causalidad”. La expresión para el coeficiente de correlación de Pearson (o r de Pearson) es:
r = \frac{\sum (x_i - \bar x)(y_i - \bar y)}{\sqrt{\sum(x_i - \bar x)^2 \sum(y_i - \bar y)^2}}

Veamos ahora algunos ejemplos en R. El comando para obtener la r de Pearson (coeficiente de correlación) es cor(x,y):

> x = seq(0, 10, by=0.1)
> y = 5*x + 7
> plot(x,y)
> cor(x,y)
[1] 1
> y = -5*x + 7
> plot(x,y)
> cor(x,y)
[1] -1

Agreguemos un poco de ruido (incertidumbre) a este ejemplo:

> x = seq(0, 10, by=0.1)
> y = 5*x + 7
> ruido = rnorm(length(x), mean=0, sd=3)
> plot(x, y + ruido)
> cor(x, y + ruido)
[1] 0.9818713

Vemos que el coeficiente de correlación (r de Pearson, en este caso) ya no es exactamente 1 debido al ruido introducido; sin embargo, aún existe una alta asociación (0.9818) entre las variables.

Los comandos importantes para regresión lineal son lm (linear model), abline. Se ejemplifican a continuación:

> x = seq(0, 10, by=0.1)
> y = 5*x + 7
> ruido = rnorm(length(x), mean=0, sd=3)
> plot(x, y + ruido)
> lm(y~x)

Call:
lm(formula = y ~ x)

Coefficients:
(Intercept)            x  
          7            5  

> lm(y + ruido ~ x)

Call:
lm(formula = y + ruido ~ x)

Coefficients:
(Intercept)            x  
      7.498        4.913  

> yr = y + ruido
> summary(lm(yr ~ x))

Call:
lm(formula = yr ~ x)

Residuals:
    Min      1Q  Median      3Q     Max 
-6.0060 -1.7720 -0.4454  1.8719  6.3629 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   7.4980     0.5383   13.93   <2e-16 ***
x             4.9134     0.0930   52.83   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 2.725 on 99 degrees of freedom
Multiple R-squared: 0.9657,     Adjusted R-squared: 0.9654 
F-statistic:  2791 on 1 and 99 DF,  p-value: < 2.2e-16 

> plot(x,yr)
> abline(lm(yr ~ x))
> plot(x,yr)
> abline(7.4980, 4.9134)

regr-lin


Prueba t de una muestra

Empecemos dejando claro que la prueba t está diseñada en el supuesto que los datos vienen de una distribución normal N(\mu, \sigma^2). Podemos entonces estimar los parámetros \mu, \sigma de la distribución muestral a través de la media \overline{x} y desviación estándar s. La prueba t se emplea para probar hipótesis acerca de las medias de dos muestras, e.g. para probar si existe alguna diferencia entre el grupo control y el grupo experimental. Veamos algunos ejemplos.

Dieta recomendada.

Unos investigadores le han pedido a 11 estudiantes de secundaria que registren sus hábitos alimenticios. En base a las observaciones de los estudiantes se tienen los siguientes datos de consumo energético promedio en kilo-joules (1 kJ = 239 Cal). La dieta recomendada tiene un valor energético de 7725 kJ. ¿Qué tanto se desvían las dietas de las estudiantes de la dieta recomendada?

> consumo = c(5260,5470,5640,6180,6390,6515,6805,7515,7515,8230,8770)
> summary(consumo)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   5260    5910    6515    6754    7515    8770 
> t.test(consumo, mu=7725)

        One Sample t-test

data:  consumo 
t = -2.8208, df = 10, p-value = 0.01814
alternative hypothesis: true mean is not equal to 7725 
95 percent confidence interval:
 5986.348 7520.925 
sample estimates:
mean of x 
 6753.636 

La primera línea que contiene información que necesitamos interpretar es:

t = -2.8208, df = 10, p-value = 0.01814

En alguna época hubiéramos tenido que calcular t y luego buscar en una tabla de una distribución t el valor de p. Acá podemos leerlo inmediatamente. Dado que p < 0.05, tenemos entonces que rechazar la hipótesis nula: la dieta promedio de los estudiantes es igual a la dieta recomendada (ojo a la línea que despliega la hipótesis alternativa) Nota: un intervalo de confianza siempre se calcula seleccionando primero un nivel de confianza (parámetro conf.level por defecto es 0.95).


Consumo de gasolina.

La compañía de carros X ha anunciado que la nueva SUV tiene un consumo de gasolina de 17 millas por galón. Los editores de una revista automovilística consideran que la compañía X ha mentido y que dicha SUV hace menos millas por galón. Toman las siguientes mediciones y corren una prueba t:

> mpg = c(11.4,13.1,14.7,14.7,15.0,15.5,15.6,15.9,16.0,16.8)
> t.test(mpg, mu=17, alt='less')

        One Sample t-test

data:  mpg 
t = -4.2847, df = 9, p-value = 0.001018
alternative hypothesis: true mean is less than 17 
95 percent confidence interval:
     -Inf 15.78127 
sample estimates:
mean of x 
    14.87 

Dado que p-value < 0.05, nuevamente rechazamos la hipótesis nula ("la nueva SUV hace 17 millas por galón") en favor de la hipótesis alternativa ("la nueva SUV hace menos de 17 millas por galón"). Noten el parámetro alt='less' que modifica el análisis y la hipótesis alternativa.


Costos de cafetería.

La cafetería de la UVG dice que, en promedio, los estudiantes gastan Q.101.75 semanales en refacciones, almuerzos, etc. Un grupo de estudiantes de estadística realiza un muestreo aleatorio y encuentran los siguientes consumos semanales: 140, 125, 150, 124, 143, 170, 125, 94, 127, 53. La hipótesis nula en este caso es que el consumo semanal por estudiante es de Q101.75. La hipótesis alternativa es que el consumo semanal es mayor a lo que dice la cafetería.

Nota: dado que la prueba t asume que los datos vienen de una distribución normal, lo primero que deberíamos hacer es visualizar qué tan normales son los datos. Pueden emplear el comando qqnorm (vean su documentación). Un ejemplo:

> qqnorm(rnorm(1000))

Ahora la corrida de la prueba t:

> gasto = c(140, 125, 150, 124, 143, 170, 125, 94, 127, 53)
> qqnorm(gasto)
> t.test(gasto, mu=101.75, alt='greater')

        One Sample t-test

data:  gasto 
t = 2.291, df = 9, p-value = 0.02385
alternative hypothesis: true mean is greater than 101.75 
95 percent confidence interval:
 106.4169      Inf 
sample estimates:
mean of x 
    125.1 

Dado el valor p obtenido (i.e., p-value < 0.05), nuevamente rechazamos la hipótesis nula; en otras palabras: en base a los datos recopilados, el consumo semanal (por estudiante) en la cafetería es mayor a Q101.75. Nuevamente noten el parámetro alt='greater' y cómo modifica la hipótesis alternativa.


Análisis de Varianza

Supongamos que tenemos varias mediciones separados en distintos grupos. Empecemos con un poco de notación: x_{i,j} representa entonces la observación (o medición) j en el grupo i. Ahora un poco de manipulación algebraica:
\displaystyle x_{i,j} = x_{i,j}
\displaystyle = \bar x_i + x_{i,j} - \bar x_i
\displaystyle = \bar x + \bar x_i - \bar x + x_{i,j} - \bar x_i
\displaystyle = \bar x + (\bar x_i - \bar x) + (x_{i,j} - \bar x_i)

donde (\bar x_i - \bar x) es la desviación de la media del grupo (\bar x_i) con respecto a la media global (\bar x); (x_{i,j} - \bar x_i) es la desviación de la observación con respecto a la media del grupo. Esto corresponde entonces al modelo estadístico:
\displaystyle X_{i,j} = \mu + \alpha_i + \varepsilon_{i,j}
donde \varepsilon_{i,j} \sim N(0, \sigma^2). Tenemos entonces las siguientes definiciones:

  • Variación dentro de grupos: \displaystyle SS_W = \sum_i \sum_j (x_{i,j} - \bar x_i)^2
  • Variación entre grupos: \displaystyle SS_B = \sum_i \sum_j (\bar x_i - \bar x)^2

La parte clave del desarrollo matemático está en notar que la varianza total puede escribirse como la suma de las varianzas anteriores:
\displaystyle SS = \sum_i \sum_j (x_{i,j} - \bar x)^2 = SS_W + SS_B
En palabras, la agrupación de los datos explica parte de la varianza, por lo que una agrupación informativa explica una gran parte de la varianza. Nota: piensen en el análisis de varianza como una generalización al análisis de prueba t.
Definimos finalmente la estadística F como:
\displaystyle F = \frac{SS_B/(k-1)}{SS_W/(N-k)}
para N observaciones separadas en k grupos.

El comando en R que empleamos para el análisis de varianza es oneway.test. Veamos un ejemplo en donde tenemos observaciones separadas en 4 grupos (datos1, datos2, datos3, datos4). En este caso las observaciones son datos simulados empleando una distribución normal con $\latex \mu = 75, \sigma = 20$:

> datos1 = rnorm(20, mean=75, sd=20)
> datos2 = rnorm(20, mean=75, sd=20)
> datos3 = rnorm(20, mean=75, sd=20)
> datos4 = rnorm(20, mean=75, sd=20)
> 
> datos.tabla = data.frame(datos1, datos2, datos3, datos4)
> boxplot(datos.tabla)
> 
> datos.stack = stack(datos.tabla)
> names(datos.stack)
[1] "values" "ind"   
> 
> oneway.test(values ~ ind, data=datos.stack, var.equal=T)

        One-way analysis of means

data:  values and ind 
F = 0.9687, num df = 3, denom df = 76, p-value = 0.412

anova

Dado que p > 0.05 aceptamos la hipótesis nula que las medias de todos los grupos son iguales; en símbolos \bar x_1 = \bar x_2 = \cdots = \bar x_k.

Nota: Otros comandos que proporcionan resultados equivalentes a oneway.test son:

  • summary(aov(values ~ ind, data=datos.stack))
  • anova(lm(values ~ ind, data=datos.stack))
Tagged with: , ,

Obteniendo resultados estadísticos con R, parte III

Posted in Tutoriales, Uncategorized by hector on abril 28, 2009

Esta es la tercera parte de una breve introducción a los comandos del sistema estadístico R. El objetivo hasta ahora ha sido obtener fácilmente resultados de estadísticas descriptivas; en esta parte nos concentraremos en desarrollar simulaciones de probabilidades.


Recordemos los problemas de la Hoja de Trabajo 5…

Disco duro I. Considere la posición inicial de un disco duro como cero (0) y la posición final como uno (1). La cabeza del disco duro lee archivos en posiciones aleatorias del disco. Luego de cada lectura, la cabeza del disco regresa a la posición inicial, i.e. cero (0). ¿Cuál es la distancia promedio que se mueve la cabeza del disco duro? Ayuda: asuma una distribución uniforme.

En clase respondimos esta pregunta de dos formas: la manera analítica y la manera computacional. Sin embargo, en ese entonces empleamos el lenguaje Python por tener una sintaxis clara y porque algunos de ustedes ya lo habían utilizado antes. Recordemos que en Python el script se veía así:

from pylab import *

def disco1(N=1000):
    """Resuelve el Problema 1 de la Hoja de Trabajo 5.

    Input: un numero grande (N) de lecturas en el disco.
    """
    dist_acum = 0   # distancia acumulada
    for k in range(N):   # loop para un numero grande (N) de lecturas
        pos = random()   # posicion de lectura
        dist = 2*pos     # distancia de ida y vuelta
        dist_acum = dist_acum + dist   # acumula la distancia
    return dist_acum / N   # normaliza para el numero de lecturas realizadas

En una terminal de Python podemos cargar el archivo con el comando run, de la siguiente manera:

In [23]: run 2009-cu109-hoja05.py

In [24]: disco1(10000)
Out[24]: 1.0039466111356798

In [25]: disco1(10000)
Out[25]: 1.0055115295663883

In [26]: disco1(10000)
Out[26]: 0.99116255989859414

In [27]: disco1(10000)
Out[27]: 1.0023458574057422

Vemos también que varias corridas de la función disco1 nos devuelven valores que tienden a 1, la respuesta esperada.

Ahora bien, en el sistema estadístico R el script no cambia mucho. Antes de ver el script completo, vale la pena mencionar un caso sencillo del enunciado de iteración for en R:

> seq(5)
[1] 1 2 3 4 5
> for (k in seq(5)) print(k*k)
[1] 1
[1] 4
[1] 9
[1] 16
[1] 25

Acá está el script en R:

disco1 = function(N)
# Resuelve el Problema 1 de la Hoja de Trabajo 5
{
dist_acum = 0   # distancia acumulada
for (k in seq(N))   # loop para un numero grande (N) de lecturas
    {
    pos = runif(1)   # posicion de lectura
    dist = 2*pos     # distancia de ida y vuelta
    dist_acum = dist_acum + dist   # acumula la distancia
    }
return(dist_acum/N)   # normaliza para el numero de lecturas realizadas
}

Pregunta 0.

Explique cómo funciona el comando runif en el código anterior.

En una terminal de R podemos cargar el archivo con el comando source (noten la similitud con el comando run de Python). Además, varias corridas de la funcion disco1 (con un número grande de lecturas) nos devuelven nuevamente valores cercanos a 1.

> source('2009-cu109-hoja05.R')
> disco1(100000)
[1] 1.000558
> disco1(100000)
[1] 1.000702
> disco1(100000)
[1] 0.9997121
> disco1(100000)
[1] 1.003375

Aclaración: ustedes obtendrán respuestas ligeramente distintas debido a la naturaleza de los algoritmos que generan números aleatorios.

Nota: asegúrense de que directorio actual sea el adecuado, i.e. File / Change dir…, o bien, en la terminal de R pueden usar los comandos getwd y
setwd (vean la documentación de éstos para ver cómo funcionan).


Disco duro II. Considere el escenario del problema anterior con la variante que la cabeza del disco duro no regresa a la posición inicial (cero) luego de cada lectura. ¿Cuál es la distancia promedio que se mueve la cabeza del disco duro ahora?

De manera similar al problema Disco duro I, respondemos acá con una comparación entre lenguajes: Python versus R.

El script en Python:

from pylab import *

def disco2(N=1000):
    """Resuelve el Problema 2 de la Hoja de Trabajo 5.

    Input: un numero grande (N) de lecturas en el disco.
    """
    dist_acum = 0   # distancia acumulada
    pos_antes = random()   # posicion anterior de lectura
                           # es igual que empezar en la posicion cero (0)
    for k in range(N):   # loop para un numero grande (N) de lecturas
        pos = random()   # posicion actual de lectura
        dist = abs(pos - pos_antes)   # distancia entre posiciones
        dist_acum = dist_acum + dist   # acumula la distancia
        pos_antes = pos   # la posicion actual es ahora la pos. anterior
    return dist_acum / N

En una terminal de Python, las respuestas de la función disco2 tienden a 0.333 (la respuesta teórica es 1/3):

In [36]: run 2009-cu109-hoja05.py

In [37]: disco2(10000)
Out[37]: 0.33279175920353116

In [38]: disco2(10000)
Out[38]: 0.33385206536325829

In [39]: disco2(10000)
Out[39]: 0.33342096518733244

Nuevamente el script en R se ve bastante similar:

disco2 = function(N)
# Resuelve el Problema 2 de la Hoja de Trabajo 5
{
dist_acum = 0   # distancia acumulada
pos_antes = runif(1)   # posicion anterior de lectura
                       # podemos tambien empezar en cero (0)
for (k in seq(N))   # loop para un numero grande (N) de lecturas
    {
    pos = runif(1)    # posicion actual de lectura
    dist = abs(pos - pos_antes)     # distancia entre posiciones
    dist_acum = dist_acum + dist    # acumula la distancia
    pos_antes = pos   # la posicion actual es ahora la pos. anterior
    }
return(dist_acum/N)   # normaliza para el numero de lecturas realizadas
}

En una terminal de R los resultados obtenidos son similares a los siguientes:

> source('2009-cu109-hoja05.R')
> disco2(100000)
[1] 0.3333376
> disco2(100000)
[1] 0.333403
> disco2(100000)
[1] 0.3327578

Pregunta 1.

Disco duro III. Considere la posición inicial de un disco duro como cero (0) y la posición final como uno (1). La cabeza del disco duro lee archivos en posiciones aleatorias del disco. Luego de cada lectura, la cabeza del disco regresa a su posición de reposo. Si dicha posición de reposo es cero (0) tenemos el mismo problema que en Disco duro I. ¿Cuál es la distancia promedio que se mueve la cabeza del disco duro dadas las siguientes posiciones de reposo? Ayuda: nuevamente asuma una distribución uniforme.

  1. Posición de reposo en 1.
  2. Posición de reposo en 0.5.
  3. Posición de reposo en 0.25.
  4. Posición de reposo en 0.75.

Ayuda: Para responder esta pregunta debe escribir una función en R en base a la siguiente función en Python:

def disco3(N, pos_reposo):
    dist_acum = 0
    for k in range(N):
        pos = random()
        dist = 2*abs(pos_reposo - pos)
        dist_acum = dist_acum + dist
    return dist_acum / N

Nota: no olviden comentar cada línea adecuadamente (como en los ejemplos anteriores).


Pregunta 2.

En la pregunta anterior, ¿encuentran algún patrón al cambiar la posición de reposo? ¿cómo cambia la distancia promedio al barrer dicha posición de reposo del inicio (0) hasta el final (1)? ¿será acaso una función lineal, una función cuadrática, o talvez una función exponencial?

Tagged with: , ,

Obteniendo resultados estadísticos con R, parte II

Posted in Tutoriales, Uncategorized by hector on abril 22, 2009

Esta es la segunda parte de una breve introducción a los comandos del sistema estadístico R. El objetivo es obtener fácilmente resultados de estadísticas descriptivas, y simulaciones de probabilidades, entre otros.

Al empezar a aprender un nuevo sistema puede resultar muy útil tener una referencia breve a la mano (un chivo, en buen chapín). Existen varias de estas para R, pero creo que ésta es muy buena.


Manejo de Variables

Empecemos hablando un poco del manejo de variables en R. Para ingresar unos pocos datos podemos usar el comando c() (concatenar). A continuación definimos una variable pesos con los pesos en libras de una familia con 5 miembros.

> pesos = c(136, 117, 215, 170, 152)
> sum(pesos)
[1] 790
> length(pesos)
[1] 5
> promedio = sum(pesos) / length(pesos)
> promedio
[1] 158
> mean(pesos)
[1] 158

Encontramos la suma de los pesos con el comando sum, la cantidad de pesos ingresados con el comando length (longitud de la variable), y encontramos el peso promedio de la familia de dos maneras: con la secuencia de comandos sum(pesos) / length(pesos) y más fácilmente con el comando mean. R ya incluye muchos comandos para que no tengamos que reinventar la rueda (vean la “reference card” al inicio de este artículo).

Algo que vale la pena entender es el manejo de índices en R. Veamos:

> pesos
[1] 136 117 215 170 152
> pesos[1]
[1] 136
> pesos[5]
[1] 152
> pesos[1:3]
[1] 136 117 215
> pesos[pesos > 150]
[1] 215 170 152
> pesos[pesos  pesos > 150
[1] FALSE FALSE  TRUE  TRUE  TRUE
> pesos < 150
[1]  TRUE  TRUE FALSE FALSE FALSE

Regresemos ahora a los datos de la encuesta de estudiantes que tratamos en la parte I, podemos clasificar (o filtrar) los datos empleando índices:

> Zap.mujer = Zapatos[Sexo == 'mujer']
> Zap.hombre = Zapatos[Sexo == 'hombre']
> boxplot(Zapatos, Zap.mujer, Zap.hombre, names=c('Todos', 'Mujeres', 'Hombres'))

zapatos

Parece cumplirse el estereotipo de que los hombres tienen pocos pares de zapatos, mientras que las mujeres tienen muchos.

Una forma más inmediata de filtrar los datos es empleando la sintaxis y~x; que se lee: la variable y descrita por la variable x.

> boxplot(Zapatos ~ Sexo)
> boxplot(Estatura ~ Sexo)

zapatos2

estatura

En esta última gráfica se aprecia que las estaturas de los hombres son, en general, mayores a las estaturas de las mujeres.

Resulta claro ahora que hay diferencias marcadas entre las características y hábitos de los estudiantes según el sexo.

Planteemos ahora otro tipo de pregunta: ¿estarán las horas de sueño de los estudiantes relacionadas con la hora a la cual deciden acostarse usualmente? Podemos ver la relación entre las variables Acostarse y horas.de.sueno (recuerden que habíamos definido horas.de.sueno = Levantarse - Acostarse):

> horas.de.sueno = Levantarse - Acostarse
> plot(Acostarse, horas.de.sueno)

scatter

Esta gráfica indica que los estudiantes que se acuestan más tarde, duermen menos.


Preguntas

  1. En promedio, ¿qué tanto más altos son los hombres que las mujeres entre los estudiantes encuestados?
  2. Elabore diagramas de caja y bigotes paralelos (boxplot) para los costos de corte de pelo, clasificados por sexo. ¿Qué puede decir de la gráfica obtenida?
  3. Construya una gráfica de barras de las frecuencias de los distintos valores que toma la variable Numero. Recuerde que esta variable describe un número aleatorio entre 1 y 10 proporcionado por el estudiante. ¿Por qué cree que, al tratar de pensar en un número aleatorio, algunos números son evitados sistemáticamente mientras otros ocurren frecuentemente?
Tagged with: , ,

Buscador para R

Posted in Uncategorized by hector on abril 21, 2009

rlogoAlgunas veces Google tiene problemas con ciertos términos de búsqueda, por ejemplo: R, latex. En estos casos es mejor usar buscadores especializados como RSeek. RSeek es un buscador en torno a R y temas estadísticos en general. Pueden usarlo para encontrar funciones, documentación, etc.

Tagged with: ,

Obteniendo resultados estadísticos con R, parte I

Posted in Tutoriales, Uncategorized by hector on abril 20, 2009

Esta es una breve introducción a los comandos del sistema estadístico R. El objetivo es obtener fácilmente resultados de estadísticas descriptivas, y simulaciones de probabilidades, entre otros.

Introducción

Una sesión interactiva de R se ve así:
screenshot1
Es acá en donde estaremos ingresando comandos e interactuando con los resultados obtenidos.

Ingreso de Datos

Lo primero que deben hacer es bajar el archivo estudiantes.doc. Este es un archivo de texto que contiene los resultados de una encuesta suministrada a 657 estudiantes. Las primeras líneas de este archivo se ven así:

"Estudiante","Estatura","Sexo","Zapatos","Numero","DVDs","Acostarse","Levantarse","CortePelo","Trabajo","Bebida"
1,67,"mujer",10,5,10,-2.5,5.5,60,30,"agua"
2,64,"mujer",20,7,5,1.5,8,0,20,"gaseosa"
3,61,"mujer",12,2,6,-1.5,7.5,48,0,"leche"
4,61,"mujer",3,6,40,2,8.5,10,0,"agua"
5,70,"hombre",4,5,6,0,9,15,17.5,"gaseosa"

Las 11 columnas anteriores corresponden a 11 variables descritas a continuación:

Estudiante: número de estudiante
Estatura: estatura en pulgadas
Sexo: sexo del estudiante
Zapatos: número de pares de zapatos que posee
Numero: un número entre 1 y 10 proporcionado por el estudiante
DVDs: número de películas en DVD que posee
Acostarse: hora a la que usualmente se acuesta
Levantarse: hora a la que usualmente se levanta
CortePelo: costo del último corte de pelo
Trabajo: número de horas que trabaja en una semana
Bebida: bebida preferida (agua, gaseosa, leche)

Antes de poder analizar estos datos debemos importarlos a una sesión de R. Uno de los comandos para importar es read.table (pueden ver la ayuda ejecutando ?read.table, ésta es en general la forma de obtener ayuda en R: ?comando). Además, asignaremos esta tabla de datos a la variable estudiantes. Todo ésto en una sola línea:

> estudiantes = read.table('estudiantes.doc', sep=',', header=TRUE)

Los argumentos en un comando nos permiten ajustar su uso. En el comando anterior el argumento sep=',' nos permite especificar que los distintos valores en el archivo están separados por comas; mientras que el argumento header=TRUE toma en cuenta el encabezado de las columnas.

Pueden ahora ver las respuestas que dió el primer estudiante (atención a la notación):

> estudiantes[1,]
Estudiante Estatura  Sexo Zapatos Numero DVDs Acostarse Levantarse CortePelo
1          1       67 mujer      10      5   10      -2.5        5.5        60
Trabajo Bebida
1      30   agua

O bien, las respuestas dadas por los estudiantes 17 al 20:

> estudiantes[17:20,]
Estudiante Estatura   Sexo Zapatos Numero DVDs Acostarse Levantarse
17         17       62  mujer       8      7   NA       0.0        8.5
18         18       66 hombre       3      7   20       2.5        8.0
19         19       71 hombre       6      7    0       0.5        7.5
20         20       64  mujer       4      7    8      -1.5        7.5
CortePelo Trabajo  Bebida
17        40     0.0    agua
18        13     0.0 gaseosa
19        12    12.5    agua
20        30    16.0    agua

Para cargar las variables de la tabla (tabla se dice data frame en el lingo de R) al ambiente se utiliza el comando attach:

> attach(estudiantes)

Ahora podemos referirnos a las variables incluidas en la data frame (tabla) estudiantes; por ejemplo, las estaturas de los primeros 10 estudiantes o el tamaño de la variable Estatura (657 datos correspondientes a cada uno de los estudiantes encuestados):

> Estatura[1:10]
[1] 67 64 61 61 70 63 61 64 66 65
> length(Estatura)
[1] 657

Descripción de Datos

El comando table permite contar rápidamente los distintos valores que toma una variable, por ejemplo:

> table(Bebida)
Bebida
agua gaseosa   leche
355     178     113
> table(Sexo)
Sexo
hombre  mujer
222    435

Noten que la mayoría de los estudiantes prefieren tomar agua (355); podríamos acaso decir que esta es una población estudiantil saludable? Será porque la mayoría son mujeres (435)?

Ahora bien, podemos graficar estos resultados con el comando barplot:

> t = table(Bebida)
> barplot(t)

barplot

En una sesión de R los comandos son interactivos, como en muchos otros sistemas de cómputo, e.g. Matlab, Python, Octave, Scilab, etc. Sin embargo, si queremos exportar los gráficos obtenidos podemos seguir una secuencia de comandos similar a la siguiente:

> png('grafica-01.png')
> barplot(t)
> dev.off()

En Windows también pueden hacerlo con el menú, como en la siguiente figura:
screenshot2

Supongan ahora que estamos interesados en las horas usuales de sueño de los estudiantes; pues consideramos que un estudiante cansado es un estudiante que rinde menos. No tenemos esta información de manera explícita, pero sí la tenemos de manera implícita, i.e. simplemente comparamos las variables Acostarse y Levantarse; definimos la variable horas.de.sueno como la diferencia:

> horas.de.sueno = Levantarse - Acostarse
> summary(horas.de.sueno)
Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's
2.500   6.500   7.500   7.385   8.500  12.500   4.000

Inmediatamente ejecutamos el comando summary que devuelve un resumen numérico de los datos, i.e. mínimo, máximo, primer cuartil, tercer cuartil, mediana (segundo cuartil), así como la media. Vemos entonces que estos estudiantes usualmente duermen alrededor de 7.5 horas.

Noten que al definir la variable horas.de.sueno como la resta Levantarse - Acostarse, R ejecuta esta operación componente por componente; es decir, de forma vectorizada.

Gráficamente podemos ver las horas de sueño de varias maneras. Una alternativa es usar un diagrama de caja y bigotes (“box-and-whisker”).
horas-de-sueno

Otra alternativa para visualizar las horas de sueño es emplear un histograma:
horas-de-sueno2

Los comandos que generaron estas gráficas son:

> boxplot(horas.de.sueno)
> hist(horas.de.sueno)

¿A que distribución de probabilidad les recuerda este histograma de horas de sueño?

Preguntas

  1. Utilice el comando summary para obtener un resumen numérico de la variable DVDs.
  2. Construya un histograma de frecuencias para DVDs.
  3. Construya una gráfica de barras de las frecuencias de los distintos valores que toma la variable DVDs. Ayuda: emplee los comandos table y barplot.
Tagged with: , ,

Examen Parcial 1, solución

Posted in Uncategorized by hector on abril 20, 2009

Acá está el Examen Parcial 1, así como la solución del mismo. Por favor estudien este material.

Nota: El último día para solicitar revisión de este examen es el lunes 27 de Abril.

Tagged with: ,

test

Posted in Uncategorized by hector on abril 20, 2009

test \int_{a}^{b}f(x)~dx

for k in [1,2,3]:
    print k
> table(t)
> hist(t)
> boxplot(t)
> table(t)
> hist(t)
> boxplot(t)
Seguir

Recibe cada nueva publicación en tu buzón de correo electrónico.