Fundamentos Estadísticos

Utilizando el Software R

Lic. Juan Isaula

Introducción a R

¿Qué es R?

  • R es un lenguaje y un ambiente para el manejo de datos, cálculos, y gráficos en código libre.

  • Además de ser gratuito, los desarrollos en R se actualzan más rápido que cualquier otro de los costosos software comerciales que se encuentran en el mercado.

Ventajas de usar R


  • Es software libre y por tanto su costo es nulo.

  • Es multiplataforma: existen versiones para LinuX, Mac y Windows. Los procedimientos y análisis

  • Capacidad para acceder a datos en múltiples formatos. Dispone de librerías para leer datos desdes SPSS, SAS, Access, MySQL, Excel, etc.

  • Generación de gráficos de alta calidad.

  • Todos los algoritmos implementados en R pueden ser vistos e interpretados por cualquier usuario, por lo que este puede saber exactamente que es lo que hace el ordenador cuando ejecuta un comando.

Desventajas de Usar R


  • Hay empresas que por políticas no pueden instalar software libre en sus maquinas cada una tiene su politica, sus software de preferencia, sus necesidades, etc.

  • Una de las principales desventajas es que hasta hace poco el uso de R estaba limitado a entornos universitarios y de usuarios con gran conocimiento de la estadística y la programación.

  • No hay nadie a quien reclamar si algo falla, ni hay un departamento de atención al cliente que nos diga qué podemos hacer si algo va mal, si alguién procedimiento nos da un error, o simplemente si no sabemos qué sintaxis utilizar.

Instalación de R

Para realizar la instalación de R y RStudio en windows, Mac, Ubuntu o Linux se debe ingresar a los siguientes sitios web:

Operaciones Básicas

Operador Matemático Función en R
\(\sqrt{x}\) sqrt(x)
\(e^x\) exp(x)
\(x!\) factorial(x)
\(\log(x)\) log(x)
\(\sin(x)\) sin(x)
\(\cos(x)\) cos(x)
\(\tan(x)\) tan(x)
\(\pi\) pi
  • <- asignar valores a una variable.
w <- 7

Vectores

Un vector es una secuencia ordenada de datos, los cuales han de ser del mismo tipo.

Tipos de datos que se pueden almacenar en un vector se destacan los siguientes:

  • logical (TRUE o FALSE)

  • integer (números enteros)

  • numeric (números reales)

  • character (palabras)

vector <- c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9)
vector
 [1] 0 1 2 3 4 5 6 7 8 9

vectorletra <- c("a", "b", "c", "d", "e")
vectorletra
[1] "a" "b" "c" "d" "e"
vectorpalabra <- c("Fundamentos", "Estadística", "en", "R")
vectorpalabra
[1] "Fundamentos" "Estadística" "en"          "R"          
w <- c(0:10)
w
 [1]  0  1  2  3  4  5  6  7  8  9 10
y <- seq(0, 100, by = 10)
y
 [1]   0  10  20  30  40  50  60  70  80  90 100

Operaciones con Vectores

w + 5
 [1]  5  6  7  8  9 10 11 12 13 14 15
10*w
 [1]   0  10  20  30  40  50  60  70  80  90 100
sqrt(w)
 [1] 0.000000 1.000000 1.414214 1.732051 2.000000 2.236068 2.449490 2.645751
 [9] 2.828427 3.000000 3.162278
w^2
 [1]   0   1   4   9  16  25  36  49  64  81 100

Otras Funciones para aplicar a vectores

  • max() y min() calculan sus valores máximos y mínimos respectivamente.

  • sum() calcula la suma

  • prod() calcula el producto

  • mean() calcula la media

  • sort() ordena los elementos del vector en el orden natural creciente del tipo de datos que lo forman, se puede incluir en su argumento el parámetro decreasing = TRUE.

w
 [1]  0  1  2  3  4  5  6  7  8  9 10

Realizando operaciones con estas nuevas funciones se tiene:

max(w)
[1] 10
min(w)
[1] 0
sum(w)
[1] 55
prod(w)
[1] 0
mean(w)
[1] 5
sort(w)
 [1]  0  1  2  3  4  5  6  7  8  9 10

Matrices

Una matriz no es más que un arreglo de números en m renglones y n columnas.

La función matrix(vector, nrow = n,) define una matriz de n filas.

matrix(1:6, nrow = 2)
     [,1] [,2] [,3]
[1,]    1    3    5
[2,]    2    4    6
matrix(1:6, nrow = 3)
     [,1] [,2]
[1,]    1    4
[2,]    2    5
[3,]    3    6
matrix(1:6, nrow = 2, byrow = TRUE)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6

Otras formas de crear Matrices

La función

rbind(vector1, vector2, vector3, . . . )

permite la construcción de la matriz de filas

rbind(c(1, 0, 2), c(2, 3, 6), c(1, 2, 0))
     [,1] [,2] [,3]
[1,]    1    0    2
[2,]    2    3    6
[3,]    1    2    0

cbind(): obtiene la matriz cuyas columnas son los vectores a los que se aplica.

     [,1] [,2] [,3]
[1,]    1    2    1
[2,]    0    3    2
[3,]    2    6    0

DataFrames

Un dataframe tiene la apariencia de una matriz, pero con la diferencia que cada columna de un dataframe tiene la apariencia de una matriz. Particularmente cada columna de un dataframe puede contener datos de un tipo diferente.

Creación de un DataFrame

Se usa la función data.frame(v1, v2, . . .) aplicada a los vectores en el orden en el que queramos disponer las columnas de la tabla.

Ejemplo

# Número  de polizas
Poliza <- c(1:9)

# Tipo de poliza 
Tipo <- c("Hogar", "Auto", "Auto", "Auto", "Hogar", "Hogar", "Auto", 
          "Auto", "Hogar")

# Fecha de creación de la póliza
Fecha <- c("12/12/2016", "08/02/2014", "10/08/2012", "01/01/2015", 
           "21/11/2011", "18/01/2016", "12/04/2005", "29/03/2007",
           "18/02/2009")

# Coberturas incluidas en la poliza
Coberturas <- c("Incendio", "Robo", "Terceros", "Robo", "Robo",
                "Incendio", "Terceros", "R. Civil", "Incendio")

# Datos de los titulares de las polizas
Nombre <- c("Carlos", "Nancy", "Pedro", "Cecilia", "Ricardo", "Sofia",
            "Armando", "Vicente", "Fernando")

# Genero, sexo, edad, estado de providencia
Sexo <- c("M", "F", "M", "F", "M", "F", "M", "M", "M")

Edad <- c(25, 35, 45, 47, 24, 43, 33, 31, 40)

Estado <-  c("Campeche", "Chiapas", "Ciudad de M ́exico", "Coahuila",
            "Durango", "Guanajuato", "Guerrero", "Hidalgo", "Jalisco")

dataframe= data.frame(Poliza, Tipo, Fecha, Coberturas, Nombre, Sexo, Edad)

dataframe
  Poliza  Tipo      Fecha Coberturas   Nombre Sexo Edad
1      1 Hogar 12/12/2016   Incendio   Carlos    M   25
2      2  Auto 08/02/2014       Robo    Nancy    F   35
3      3  Auto 10/08/2012   Terceros    Pedro    M   45
4      4  Auto 01/01/2015       Robo  Cecilia    F   47
5      5 Hogar 21/11/2011       Robo  Ricardo    M   24
6      6 Hogar 18/01/2016   Incendio    Sofia    F   43
7      7  Auto 12/04/2005   Terceros  Armando    M   33
8      8  Auto 29/03/2007   R. Civil  Vicente    M   31
9      9 Hogar 18/02/2009   Incendio Fernando    M   40

Instalación de Librerías

  • installed.packages("name") \(\Rightarrow\) Instalador paquetes

  • library(name) \(\Rightarrow\) Llamar o cargar el paquete

install.packages("tidyverse")

library(tidyverse)

Exportar e Importar Datos

Primero cargamos la librería tidyverse

library(tidyverse)
  • Importar datos tipo .csv

    read_csv("ruta archivo.csv")
  • Exportar datos como .csv

    write_csv("ruta", "nombre_achivo.csv")

Observacion:

tidyverse no es la única librería que nos permite exportar e importar conjuntos de datos.

  • library(xlsx) \(\Rightarrow\) R a Excel
  • library(haven) \(\Rightarrow\) R a SPSS
  • library(foreign) \(\Rightarrow\) R a STATA

Bibliografías Sugeridas

Estadística con R

Conjunto de Datos a utilizar

  • Datos sobre los hábitos de sueño de diferentes mamíferos.

  • Estan contenidos por defecto en R en el paquete de ggplot2.

library(ggplot2)
data("msleep")
head(msleep)
# A tibble: 6 × 11
  name    genus vore  order conservation sleep_total sleep_rem sleep_cycle awake
  <chr>   <chr> <chr> <chr> <chr>              <dbl>     <dbl>       <dbl> <dbl>
1 Cheetah Acin… carni Carn… lc                  12.1      NA        NA      11.9
2 Owl mo… Aotus omni  Prim… <NA>                17         1.8      NA       7  
3 Mounta… Aplo… herbi Rode… nt                  14.4       2.4      NA       9.6
4 Greate… Blar… omni  Sori… lc                  14.9       2.3       0.133   9.1
5 Cow     Bos   herbi Arti… domesticated         4         0.7       0.667  20  
6 Three-… Brad… herbi Pilo… <NA>                14.4       2.2       0.767   9.6
# ℹ 2 more variables: brainwt <dbl>, bodywt <dbl>

Descripción del DataSet

Columna Descripción
name nombre común
genus rango taxonómico
vore carnívoro, omnívoro o herbívoro
order Estado de conservación del mamífero
sleep_total cantidad total de sueño, en horas
sleep_rem sueño rem, en horas
sleep_cycle duración del ciclo de sueño
awake cantidad de tiempo que pasa despierto en horas
brainwt peso del cerebro en kilogramos
bodywt peso corporal en kilogramos

Histograma

Antes de sumergirnos en esto, recordemos cómo funcionan los histogramas.

msleep %>% 
  ggplot(aes(sleep_total)) + 
  geom_histogram() + 
  labs(title = "Distribucción del tiempo de sueño de varios mamíferos", 
       x = "Sleep_total")

Medidas de Tendencia Central

Estadística Función en R
Media mean()
Mediana median()
Moda mod()

Calculemos estas estadísticas en R para el número total de horas de sueño de los mamíferos de nuestro dataset.

# Media 
mean(msleep$sleep_total)
[1] 10.43373
# Mediana
median(msleep$sleep_total)
[1] 10.1

# Moda
msleep %>% count(sleep_total, sort = TRUE)
# A tibble: 65 × 2
   sleep_total     n
         <dbl> <int>
 1        12.5     4
 2        10.1     3
 3         5.3     2
 4         6.3     2
 5         8.4     2
 6         8.7     2
 7         9.1     2
 8         9.4     2
 9         9.8     2
10        10.3     2
# ℹ 55 more rows

Note que existen 4 mamíferos en nuestro dataset que duermen durante 12.5 horas.

¿Cómo sabemos que estadística o métrica usar?

Para responder la interrogante, consideremos el caso que tenemos todos los insectivoros de nuestra data:

library(tidyverse)
msleep %>% 
  filter(vore == "insecti") %>% 
  summarize(mean_sleep = mean(sleep_total),
            median_sleep = median(sleep_total))
# A tibble: 1 × 2
  mean_sleep median_sleep
       <dbl>        <dbl>
1       14.9         18.1

Considere ahora que hemos descubierto un nuevo insectívoro que misteriosamente nunca duerme. Con este supuesto suponga que nuestro resultados son ahora

\[ media = 11.5 \hspace{1cm}y\hspace{1cm} mediana = 17.6 \]

La media se puede reducir en más de 3 horas, mientras que la mediana cambió en menos de una hora. Esto se debe a que la media es mucho más sensible a los valores extremos que la mediana. Esto en datos simétricos

¿Sesgo?

Si los datos estan sesgados, lo que significa que no son simétricos como se muestran a continuación.

Cuando los datos están sesgados, la media y la mediana son diferentes. La media se tira en dirección del sesgo, por lo que es más baja que la mediana en los datos sesgados a la izquierda y más alta que la media en los datos sesgados a la derecha.

Ejercicio

Dataset food_consumption indice de carbono alimentario 2018 en kilogramos de alimentos consumidor por persona por año en cada país en cada categoría de alimento (consumption).

  • co2_emissions información sobre la huella de carbino de alguna categoría de alimento consumption medida en kilogramos de dióxido de carbono \(CO_2\), por persona en cada país.

  • Encontremos medidas de tendencia central para comparar el consumo de alimentos en EE.UU y Bélgica.

# Cargamos el conjunto de datos en R
food_consumption <- readRDS("food_consumption.rds")

head(food_consumption)
# A tibble: 6 × 4
  country   food_category consumption co2_emission
  <chr>     <fct>               <dbl>        <dbl>
1 Argentina pork                10.5         37.2 
2 Argentina poultry             38.7         41.5 
3 Argentina beef                55.5       1712   
4 Argentina lamb_goat            1.56        54.6 
5 Argentina fish                 4.36         6.96
6 Argentina eggs                11.4         10.5 
# Creemos dos marcos de datos: 

# 1. uno que contenga las filas de consumo de Bélgica
consumo_belgica <- food_consumption %>% 
  filter(country == "Belgium")

head(consumo_belgica, 5)
# A tibble: 5 × 4
  country food_category consumption co2_emission
  <chr>   <fct>               <dbl>        <dbl>
1 Belgium pork                38.6         137. 
2 Belgium poultry             12.2          13.1
3 Belgium beef                15.6         482. 
4 Belgium lamb_goat            1.32         46.2
5 Belgium fish                19.0          30.3

# 2. Otro que contenga las filas de consumo para USA
consumo_usa <- food_consumption %>% 
  filter(country == "USA")

consumo_usa
# A tibble: 11 × 4
   country food_category consumption co2_emission
   <chr>   <fct>               <dbl>        <dbl>
 1 USA     pork                27.6         97.8 
 2 USA     poultry             50.0         53.7 
 3 USA     beef                36.2       1118.  
 4 USA     lamb_goat            0.43        15.1 
 5 USA     fish                12.4         19.7 
 6 USA     eggs                14.6         13.4 
 7 USA     dairy              255.         363.  
 8 USA     wheat               80.4         15.3 
 9 USA     rice                 6.88         8.8 
10 USA     soybeans             0.04         0.02
11 USA     nuts                 7.86        13.9 

Cálculo de la media y mediana de kilogramos de alimentos consumidos por persona cada año año para ambos países.

  • Bélgica
mean(consumo_belgica$consumption)
[1] 42.13273
median(consumo_belgica$consumption)
[1] 12.59
  • Usa
mean(consumo_usa$consumption)
[1] 44.65
median(consumo_usa$consumption)
[1] 14.58

  • Optimizando el código previo en un sólo bloque de código.
food_consumption %>% 
  # Filtramos para Bélgica y EE.UU
  filter(country %in% c("Belgium", "USA")) %>% 
  # Agrupamos por país
  group_by(country) %>% 
  # Obtenemos resumen de consumo meduo y mediano
  summarise(mean_consumption   = mean(consumption),
            median_consumption = median(consumption))
# A tibble: 2 × 3
  country mean_consumption median_consumption
  <chr>              <dbl>              <dbl>
1 Belgium             42.1               12.6
2 USA                 44.6               14.6

  • Histograma de co2_emission para el arroz (rice).
food_consumption %>% 
  filter(food_category == "rice") %>% 
  ggplot(aes(co2_emission)) + 
  geom_histogram()

Respondamos las siguientes preguntas:

  1. ¿Observe el histograma para el arroz (rice) que acabamos de calcular, cuál de los siguientes términos describe mejor los datos?

    • Sin sesgo

    • Sesgado a la izquierda

    • Sesgado a la derecha

  2. ¿Qué medida de tendencia central resume mejor los kilogramos de \(CO_2\) para el arroz?

    • Media

    • Mediana

    • Ambos la media y la mediana

Medidas de Propagación

  • La varianza

    Mide la distancia promedio de cada punto de datos hasta la media de los datos.

Matemáticamente esto es

\[ s^2 = \frac{1}{n-1}\sum_{i=1}^n (y_i - \bar{y})^2 \Longrightarrow \mbox{Varianza muestral} \]

  • Desviación Estándar

    Es una medida que se utiliza para cuantificar la variación o dispersión de un conjunto de datos.

    \[ \sigma = \sqrt{s^2} \Longrightarrow \mbox{Desviación estándar de la muestra} \]

    • Un \(\sigma\) bajo indica que la mayor parte de los datos de una muestra tienden a estar agrupados cerca de su media.

    • Si \(\sigma\) es alto indica que los datos se extiende sobre un rango de valores más amplio.

  • Cuartiles

    Divide los datos en cuatro partes iguales

Sintaxys en R

Estadística o Métrica Función en R
\(s^2\) var()
\(\sigma\) sd()
cuartiles quantile()

Veamos como hacer esto en R.

# Varianza
var(msleep$sleep_total)
[1] 19.80568
# Desviación Estandar
sd(msleep$sleep_total)
[1] 4.450357
quantile(msleep$sleep_total)
   0%   25%   50%   75%  100% 
 1.90  7.85 10.10 13.75 19.90 

Gráfico de cajas para cuartiles 

ggplot(msleep, aes(y = sleep_total))+
  geom_boxplot()

Probabilidad y Números Aleatorios

¿Qué son las probabilidades?

\[P(Evento) = \frac{\mbox{# de formas en que puede ocurrir el evento}}{\mbox{# total de resultados posibles}}\]

  • La probabilidad siempre está entre cero y 100%.

  • ¿Probabilidad de que Brian sea seleccionado?

\[ P(Brian) = \frac{1}{4} = \fbox{25%} \]

Muestreo desde un marco de datos

Podemos recrear el escenario anterior de la selección de las cuatro personas, utilizando la función sample_n() de dplyr.

  • sample_n() toma de un marco de datos la cantidad de filas que queremos extrer de manera aleatoria.
vendedores <- data.frame(name = c("Amir", "Brian", "claire", "Damian"), 
                         n_ventas = c(178, 126, 75, 69))
vendedores %>% 
  sample_n(1)
    name n_ventas
1 Damian       69
vendedores %>% 
  sample_n(1)
    name n_ventas
1 claire       75

Muestreo sin reemplazo

\[ P(Brian) = \frac{1}{4} = \fbox{25%} \]

\[ P(Claire) = \frac{1}{3} = \fbox{33%} \]

Muestreo sin reemplazo\(\Rightarrow\)Dependencia

Recreemos nuestros ejemplos en R,

set.seed(5)      # valor semilla 
vendedores %>% 
  sample_n(2)
    name n_ventas
1  Brian      126
2 claire       75

Muestreo Con reemplazo\(\Rightarrow\)Independencia

Probabilidades en R

# Datos de venta de Brian
ventas_Brian <- data.frame(product = c("Product A", "Product B", "Product C", 
                                       "Product D", "Product E", "Product F", 
                                       "Product G", "Product H", "Product I",
                                       "Product J", "Product N"), 
                           n = c(23, 62, 15, 40, 5, 11, 2, 8, 7, 2, 3))

# Creemos una nueva columna llamada prob

ventas_Brian %>% 
  mutate(prob = n/sum(n))
     product  n       prob
1  Product A 23 0.12921348
2  Product B 62 0.34831461
3  Product C 15 0.08426966
4  Product D 40 0.22471910
5  Product E  5 0.02808989
6  Product F 11 0.06179775
7  Product G  2 0.01123596
8  Product H  8 0.04494382
9  Product I  7 0.03932584
10 Product J  2 0.01123596
11 Product N  3 0.01685393

Distribuciones Discretas y Continuas

Distribuciones Discretas

Se utilizan para modelar variables contables, como el lanzamiento de un dado de seis caras.

Distribuciones de probabilidad

Describe la probabilidad de cada resultado posible en un escenario

Valor Esperado

Es la media de una distribución, esto es

\[E(Y) = \sum_{y}yp(y)\]

Muestreo a partir de distribuciones discretas

# Representación de un dado justo
dado <- data.frame(n = c(1, 2, 3, 4, 5, 6))

mean(dado$n)
[1] 3.5

Tomemos una muestra de 10 lanzamientos

set.seed(4)
muestra_10 <- dado %>% 
  sample_n(10, replace = TRUE)

muestra_10
   n
1  3
2  3
3  3
4  4
5  3
6  6
7  5
8  2
9  3
10 6

Veamos la media de la muestra

mean(muestra_10$n)
[1] 3.8

Una muestra más grande

set.seed(1)
muestra_100 <- dado %>% 
  sample_n(100, replace = TRUE)

mean(muestra_100$n)
[1] 3.52

Visualizamos la muestra

ggplot(muestra_100, aes(n)) + 
  geom_histogram(bins = 6)

Ley de los Grandes Números

Lo anterior se le conoce como ley de los grandes números, ya que a medida que aumenta el tamaño de la muestra la media de la muestra se acerca a la media teórica.

Distribución Binomial

Describe la probabilidad del número de éxitos en una secuencia de ensayos independientes.

Se dice que una variable aleatoria \(Y\) tiene distribución binomial basada en \(n\) pruebas con probabilidad de éxito p si y solo si

\[ p(Y = y) =  \binom{n}{y}p^yq^{n-y}, \hspace{1.2cm} y = 0, 1,2, . . . n \hspace{0.5cm} y \hspace{0.5cm} 0\leq p \leq 1 \]

con

\[ E[Y] = n\times p \]

\[ V[Y] = np(1- p) \]

Distribución Binomial - Funciones en R

Probabilidades o simulaciones Función en R
valores aleatorios de probabilidad de una muestra dada rbinom(n, muestra, p)
\(P(Y \leq y)\) dbinom(n, muestra, p)
\(P(Y> y)\) pbinom(n, muestra, p, lower.tail = F)

Distribución Normal

Es una de las distribuciones de probabilidad más importantes que aprenderá, ya que innumerables métodos estadísticos se basan en ella y se aplica a más situaciones del mundo real que otras distribuciones.

  • Se dice que una variable aleatoria \(Y\) tiene una distribución normal de probabilidad si y sólo si para \(\sigma > 0\) y \(-\infty < \mu < \infty\),

    \[f(y) = \frac{1}{\sigma\sqrt{2\pi}}e^{-(y-\mu)^2/(2\sigma^2)}, \hspace{1cm}-\infty < y < \infty\]

con

\[ E(Y) = \mu \hspace{1cm} y \hspace{1cm} V(Y) = \sigma^2 \]

Gráficamente

Distribución Normal - Funciones en R

  • \(P(Y \leq y)\): pnorm(y, mean, sd)

  • \(P(Y > y)\): pnorm(y, mean, sd, lower.tail = F)

  • \(P(y_1 \leq Y \leq y_2)\): pnorm(y_2, mean, sd)-pnorm(y_1, mean, sd)

Caso de Estudio

Considere la data de la Encuesta Nacional de Examen de Salud y Nutrición (NHANES)

library(tidyverse)
NHANES <- read_csv("https://raw.githubusercontent.com/GTPB/PSLS20/master/data/NHANES.csv")
  • Realicemos un histograma de las alturas mayores o iguales a 135 cm de las mujeres que participaron en la encuesta.

NHANES %>% 
  filter(Gender == "female" & Height >= 135) %>% 
  ggplot(aes(x = Height)) +
  geom_histogram(aes(y = ..density..), color = "black", fill = "orange") + 
  geom_density(fill = "black", alpha = 0.2)

Tests de Normalidad

Paquetes

library(nortest) # Realiza 5 pruebas de normalidad 
library(moments) # Realiza 1 prueba de normalidad

Hipótesis

\(H_0:\) La muestra proviene de una distribución normal

\(H_1:\) La muestra no proviene de una distribución normal

Nivel de Significancia

\(\alpha = 0.05\)

  • Si \(p < \alpha\) se rechaza \(H_0\)

  • Si \(p \geq \alpha\) no se rechaza \(H_0\)

¿Nuestros datos provienen de una distribución normal?

# Pruba de Shapiro - Wilk
ad.test(NHANES$Height)

    Anderson-Darling normality test

data:  NHANES$Height
A = 452.86, p-value < 2.2e-16
# Prueba de Cramer-Von Mises
# Es útil para pequeñas muestras 
x <- rnorm(100, 150, 5)
shapiro.test(x)

    Shapiro-Wilk normality test

data:  x
W = 0.97035, p-value = 0.0235

¿Cuál es la altura media y la desviación estándar?

NHANES %>% 
  summarise(media = mean(Height, na.rm = TRUE), 
            sd = sd(Height, na.rm = TRUE))
# A tibble: 1 × 2
  media    sd
  <dbl> <dbl>
1  162.  20.2

¿Qué porcentaje de mujeres mide menos de 154 cm?

pnorm(154, mean = 161, sd = 20)
[1] 0.3631693

¿Qué porcentaje de mujeres mide más de 154 cm?

pnorm(154, mean = 161, sd = 20, lower.tail = FALSE)
[1] 0.6368307

¿Qué porcentaje de mujeres mide entre 154 - 157 cm?

pnorm(157, mean = 161, sd = 20) - pnorm(154, mean = 161, sd = 20)
[1] 0.05757094

Correlación

  • La correlación muestra la relación entre las variables. Se puede visualizar utilizando diagramas de dispersión.

    • Eje x: Variable explicativa o independiente

    • Eje y: Variable de respuesta o dependiente

    # AAPL: conjunto de datos de acciones de APPLE
    AAPL <- read_csv("AAPL.csv")
    head(AAPL)
    # A tibble: 6 × 7
      Date        Open  High   Low Close `Adj Close`     Volume
      <date>     <dbl> <dbl> <dbl> <dbl>       <dbl>      <dbl>
    1 2009-01-02  3.07  3.25  3.04  3.24        2.75  746015200
    2 2009-01-05  3.33  3.44  3.31  3.38        2.87 1181608400
    3 2009-01-06  3.43  3.47  3.30  3.32        2.82 1289310400
    4 2009-01-07  3.28  3.30  3.22  3.25        2.76  753048800
    5 2009-01-08  3.23  3.33  3.22  3.31        2.81  673500800
    6 2009-01-09  3.33  3.34  3.22  3.24        2.75  546845600

cor(AAPL$Open, AAPL$Close)
[1] 0.999731

Ambas variables están correlacionadas fuertemente.

cor(msleep$bodywt, msleep$awake)
[1] 0.3119801

La correlación entre el peso corporal y el tiempo de vigilia es solo del 0.3 , que es una relación lineal débil.

Regresión Lineal

Los Datos

# 
datos <- read_table('http://verso.mat.uam.es/~joser.berrendero/datos/EdadPesoGrasas.txt')
head(datos)
# A tibble: 6 × 3
   peso  edad grasas
  <dbl> <dbl>  <dbl>
1    84    46    354
2    73    20    190
3    65    52    405
4    70    30    263
5    76    57    451
6    69    25    302

Explorando relaciones existentes entre cada par de variables

pairs(datos)
cor(datos)
            peso      edad    grasas
peso   1.0000000 0.2400133 0.2652935
edad   0.2400133 1.0000000 0.8373534
grasas 0.2652935 0.8373534 1.0000000

Representación de la recta de mínimos cuadrados

regresion <- lm(grasas ~ edad, data = datos)
summary(regresion)

Call:
lm(formula = grasas ~ edad, data = datos)

Residuals:
    Min      1Q  Median      3Q     Max 
-63.478 -26.816  -3.854  28.315  90.881 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 102.5751    29.6376   3.461  0.00212 ** 
edad          5.3207     0.7243   7.346 1.79e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 43.46 on 23 degrees of freedom
Multiple R-squared:  0.7012,    Adjusted R-squared:  0.6882 
F-statistic: 53.96 on 1 and 23 DF,  p-value: 1.794e-07

\[ grasa = 102.58 + 5.32edad \]

Plot

ggplot(datos, aes(edad, grasas)) + 
  geom_point() +
  geom_smooth(method = "lm", se = FALSE)

El \(R^2\) o Multiple R-Squared: 0.701 mide la bondad del ajuste de la recta a los datos.

Pronósticos

Suponga que quiere predecir la cantidad de grasa para individuos de edades 30, 31, 32, . . . , 50.

nuevas_edades <- data.frame(edad = seq(30,50))
predict(regresion, nuevas_edades)
       1        2        3        4        5        6        7        8 
262.1954 267.5161 272.8368 278.1575 283.4781 288.7988 294.1195 299.4402 
       9       10       11       12       13       14       15       16 
304.7608 310.0815 315.4022 320.7229 326.0435 331.3642 336.6849 342.0056 
      17       18       19       20       21 
347.3263 352.6469 357.9676 363.2883 368.6090 
  • Para un individuo de 30 años, predecimos una cantidad de grasa de 262.2

Modelando datos de marketing

Usaremos el conjunto de datos marketing [paquete darium] Contiene el impacto de tres medios publicitarios (youtube, facebook y paper) en las ventas. Los datos son el presupuesto publicitario en miles de dólares junto con las ventas.

library(datarium)
data("marketing")
head(marketing, 4)
  youtube facebook newspaper sales
1  276.12    45.36     83.04 26.52
2   53.40    47.16     54.12 12.48
3   20.64    55.08     83.16 11.16
4  181.80    49.56     70.20 22.20

Queremos predicir las ventas futuras sobre la base del presupuesto publicitario gastado en YouTube.

pairs(marketing)
cor(marketing)
             youtube   facebook  newspaper     sales
youtube   1.00000000 0.05480866 0.05664787 0.7822244
facebook  0.05480866 1.00000000 0.35410375 0.5762226
newspaper 0.05664787 0.35410375 1.00000000 0.2282990
sales     0.78222442 0.57622257 0.22829903 1.0000000

Representación de la recta de mínimos cuadrados

ggplot(marketing, aes(youtube, sales)) +
  geom_point() +
  stat_smooth(method = lm)+
  theme_minimal()

Evaluación del Modelo

model <- lm(sales ~ youtube, data = marketing)
summary(model)

Call:
lm(formula = sales ~ youtube, data = marketing)

Residuals:
     Min       1Q   Median       3Q      Max 
-10.0632  -2.3454  -0.2295   2.4805   8.6548 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 8.439112   0.549412   15.36   <2e-16 ***
youtube     0.047537   0.002691   17.67   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3.91 on 198 degrees of freedom
Multiple R-squared:  0.6119,    Adjusted R-squared:  0.6099 
F-statistic: 312.1 on 1 and 198 DF,  p-value: < 2.2e-16

\[ ventas = 8.439112 + 0.047537\cdot youtube \]

Pronósticos

presupuesto <- data.frame(youtube = c(181.8, 53.4 ))
predict(model, presupuesto)
       1        2 
17.08127 10.97757 

Series Temporales

Definición

Una serie de tiempo es una colección de observaciones de una variable tomadas de forma secuencial y ordenada en el tiempo.

\[ y_t = \{y_1, y_2, . . . . , y_T\} \]

El análisis de series de tiempo se utiliza para muchas aplicaciones:

  • Previsión Económica

  • Pronóstico de ventas

  • Análisis presupuestario

  • Análisis del mercado de valores

  • Proyecciones de rendimiento

  • Proceso y control de calidad

  • Estudios de inventario

Ejemplo

El objetivo es el análisis, modelado y predicción de la serie temporal \(CO_2\)

# Cargamos los datos co2
co2<-read.csv("co2.csv", header = TRUE, sep = ";")

# Transformamos los datos en una serie temporal 
co2ts<-ts(co2$CO2, start = c(1959,1), frequency = 12)
print(co2ts)
        Jan    Feb    Mar    Apr    May    Jun    Jul    Aug    Sep    Oct
1959 333.13 332.09 331.10 329.14 327.36 327.29 328.23 329.55 330.62 331.40
1960 333.92 333.43 331.85 330.01 328.51 328.41 329.25 330.97 331.60 332.60
1961 334.68 334.17 332.96 330.80 328.98 328.57 330.20 331.58 332.67 333.17
1962 336.82 336.12 334.81 332.56 331.30 331.22 332.37 333.49 334.71 335.23
1963 337.95 338.00 336.37 334.47 332.46 332.29 333.76 334.80 336.00 336.63
1964 339.05 339.27 337.64 335.68 333.77 334.09 335.29 336.76 337.77 338.26
1965 341.47 341.31 339.41 337.74 336.07 336.07 337.22 338.38 339.32 340.41
1966 343.02 342.54 340.88 338.75 337.05 337.13 338.45 339.85 340.90 341.70
1967 344.28 343.42 342.02 339.97 337.84 338.00 339.20 340.63 341.41 342.68
1968 345.92 345.40 344.16 342.11 340.11 340.15 341.38 343.02 343.87 344.59
1969 347.38 346.78 344.96 342.71 340.86 341.13 342.84 344.32 344.88 345.62
1970 348.53 347.87 346.00 343.86 342.55 342.57 344.11 345.49 346.04 346.70
1971 349.93 349.26 347.44 345.55 344.21 343.67 345.09 346.27 347.33 347.82
1972 351.71 350.94 349.10 346.77 345.73                                   
        Nov    Dec
1959 331.87 333.18
1960 333.57 334.72
1961 334.86 336.07
1962 336.54 337.79
1963 337.93 338.95
1964 340.10 340.88
1965 341.69 342.51
1966 342.70 343.65
1967 343.04 345.27
1968 345.11 347.07
1969 347.23 347.62
1970 347.38 349.38
1971 349.29 350.91
1972              

#Trazamos la serie de tiempo co2ts
library(ggfortify)
autoplot(co2ts, ts.colour = "blue", ts.linetype = "dashed")
  • La serie presenta tendencia ascendente (no es estacionaria en media)

  • Es estacionaria en cuanto a la varianza, ya que no se aprecia gran variabilidad