Train Analisis

Para este ejercicio analizaremos una base de datos de trenes para analizar las tendencias, para ellos utilizaremos readr, dplyr, ggplot2 y MASS.

glimpse(bike)

Observations: 10,886
Variables: 12
$ datetime   <dttm> 2011-01-01 00:00:00, 2011-01-01 01:00:00, 2011-01-01 02:00:00, 2011-01-01 …
$ season     <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
$ holiday    <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ workingday <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ weather    <dbl> 1, 1, 1, 1, 1, 2, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 2, 2, 2, 2, 2, …
$ temp       <dbl> 9.84, 9.02, 9.02, 9.84, 9.84, 9.84, 9.02, 8.20, 9.84, 13.12, 15.58, 14.76, …
$ atemp      <dbl> 14.395, 13.635, 13.635, 14.395, 14.395, 12.880, 13.635, 12.880, 14.395, 17.…
$ humidity   <dbl> 81, 80, 80, 75, 75, 75, 80, 86, 75, 76, 76, 81, 77, 72, 72, 77, 82, 82, 88,…
$ windspeed  <dbl> 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 6.0032, 0.0000, 0.0000, 0.0000, 0.0…
$ casual     <dbl> 3, 8, 5, 3, 0, 0, 2, 1, 1, 8, 12, 26, 29, 47, 35, 40, 41, 15, 9, 6, 11, 3, …
$ registered <dbl> 13, 32, 27, 10, 1, 1, 0, 2, 7, 6, 24, 30, 55, 47, 71, 70, 52, 52, 26, 31, 2…
$ count      <dbl> 16, 40, 32, 13, 1, 1, 2, 3, 8, 14, 36, 56, 84, 94, 106, 110, 93, 67, 35, 37…

Con el comando glimpse analizamos la estructura del dataframe del archivo. Ahora filtraremos las variables que necesitamos, que son temporada, temperatura y humedad:

bike_data<-
     bike %>%  
     dplyr::select(season, atemp, humidity)
head(bike_data)

bike_data %>% 
     ggplot(aes(x=atemp,y=..density..)) + 
     geom_histogram(bins = 60) +
     geom_density() +
     facet_wrap(.~season)

Al generar estas graficas vemos la relacion entre densidad y temperatura, por estacion, en donde la estacion dos presenta dos montanas indicando que hay un vaye cuando la temperatura esta entre 20 y 25 grados, no obstante se dispara cuando la temperatura pasa los 30 grados. Asi mismo podemos observar diferentes tendencias en los datos para cada estacion, siendo la estacion 3 (otonio) la mas estable.

bike_data %>% 
     ggplot(aes(x=humidity,y=..density..)) + 
     geom_histogram(bins = 30) +
     geom_density() +
     facet_wrap(.~season)

Al analizar la humedad vemos que los datos son mas consistentes, aumentando el uso de las bicicletas cuando la humedad esta en 50 y 75, decreciendo su demananda a menos de 25 y llegando a 100.

bike_data$group<-
cut(bike_data$atemp,breaks = quantile(bike_data$atemp,probs = c(0,.25,0.5,0.75,1))) 
bike_data %>% 
     mutate(group=cut(bike_data$atemp,breaks = seq(0,30,by=0.1))) %>% 
     group_by(group) %>% 
     summarise(cnt=n()) %>% 
     ungroup() %>% 
     mutate(sum=sum(cnt)) %>% 
     mutate(prob=cnt/sum) %>% 
     dplyr::select(-sum)

Factor `group` contains implicit NA, consider using `forcats::fct_explicit_na`

bike_data

Ahora agroupamos la humedad en quartiles para que sea mas facil la visualizacion de los datos.

bike_data %>% 
     mutate(group=cut(bike_data$atemp,breaks = seq(0,30,by = 1))) %>% 
     group_by(group) %>% 
     mutate(sum(atemp))

Factor `group` contains implicit NA, consider using `forcats::fct_explicit_na`Factor `group` contains implicit NA, consider using `forcats::fct_explicit_na`Factor `group` contains implicit NA, consider using `forcats::fct_explicit_na`

bike_data

Y la secuencia de la temperatura tambien.

bike_data_season <- bike_data %>% filter(season==1)
bike_density <- kde2d(bike_data_season$atemp,bike_data_season$humidity, n=1000)
str(bike_density)

List of 3
 $ x: num [1:1000] 0.76 0.792 0.824 0.856 0.887 ...
 $ y: num [1:1000] 0 0.1 0.2 0.3 0.4 ...
 $ z: num [1:1000, 1:1000] 5.94e-15 6.37e-15 6.82e-15 7.31e-15 7.82e-15 ...

Utilizando la funcion kde2d para estimar la densidad de los datos. Ahora necesitamos convertir la densidad en una probabilidad para ello calculamos la sumatoria y luego dividimos cada valor dentro de la sumatoria, el resultado debe sumar 1.

total_density <- sum(bike_density$z)
total_density

[1] 306.9665

bike_density$z <- bike_density$z/total_density
sum(bike_density$z)

[1] 1

contour(bike_density)
text(12.2,92,"1",cex=1.5)

points(12.2,89,col="red",pch=18)
text(10.4,47,"2",cex=1.5)

points(10.4,43,col="red",pch=18)
text(21,45,"3",cex=1.5)

points(21,40,col="red",pch=18)
text(21.6,82,"4",cex=1.5)

points(21.6,78,col="red",pch=18)

hm_col_scale<-
     colorRampPalette(c("black","blue","green","orange","red"))(1000)

Ahora generamos una grafica de densidad que nos permita ver las mayores concentraciones por temporada, siendo los numeros de la grafica, maximos potenciales para cada estacion, segun la temperatura.

image(bike_density$z,  
      col = hm_col_scale,
      zlim=c(min(bike_density$z), max(bike_density$z)))
text(0.31,0.46,"Most Frequent",col="blue",cex=0.8)

points(0.31,0.45,col="blue",pch=18)

Si lo vemos como una grafica de calos vemos la mayor concentracion de probablidades segun la temperatura, la cual se encuentra entre 0.2 y 0.4.

atemp_res<-
     which(bike_density$x>10)
prob_square <- function(a,b,c,d,ked){
     atemp_res<-
          which(ked$x>a&ked$x<b)
     hum_res<-
          which(ked$y>c&ked$y<d)
     return(sum(ked$z[atemp_res,hum_res]) )
}
prob_square(10,50,50,100,bike_density)

[1] 0.4678024

Ahora utilizando una funcion cuadratica, calculamos la probabilidad de encontrar la densidad dentro de una cuadrado.

hum_res<-
     which(bike_density$y>=80&bike_density$y<=100)
hum_res

  [1]  801  802  803  804  805  806  807  808  809  810  811  812  813  814  815  816  817  818
 [19]  819  820  821  822  823  824  825  826  827  828  829  830  831  832  833  834  835  836
 [37]  837  838  839  840  841  842  843  844  845  846  847  848  849  850  851  852  853  854
 [55]  855  856  857  858  859  860  861  862  863  864  865  866  867  868  869  870  871  872
 [73]  873  874  875  876  877  878  879  880  881  882  883  884  885  886  887  888  889  890
 [91]  891  892  893  894  895  896  897  898  899  900  901  902  903  904  905  906  907  908
[109]  909  910  911  912  913  914  915  916  917  918  919  920  921  922  923  924  925  926
[127]  927  928  929  930  931  932  933  934  935  936  937  938  939  940  941  942  943  944
[145]  945  946  947  948  949  950  951  952  953  954  955  956  957  958  959  960  961  962
[163]  963  964  965  966  967  968  969  970  971  972  973  974  975  976  977  978  979  980
[181]  981  982  983  984  985  986  987  988  989  990  991  992  993  994  995  996  997  998
[199]  999 1000

sum(bike_density$z[atemp_res,hum_res])

[1] 0.1319706

hm_col_scale<-
     colorRampPalette(c("black","red"))(2)
M<-
matrix(0,ncol = 1000,nrow = 1000)
for(i in 1:1000){
     for(j in 1:1000){
           M[i,j]<-(i-250)^2/100^2+(j-250)^2/200^2 <= 1
     }
}
image(M,
      col = hm_col_scale,
      zlim=c(0, 1))

sum(bike_density$z*M)

[1] 0.06684618

Ahora con las probabilidades procedemos a dibujar un circulo. Y finalmente dibujamos un triangulo con los puntos.

x1 <- 100
y1 <- 100
x2 <- 40
y2 <- 40
x3 <- 0
y3 <- 500
m1 <- (y2-y1)/(x2 - x1)
b1 <- y1 - m1*(x1)
m2 <- (y2 - y3)/(x2 - x3)
b2 <- y2 - m2*x2
m3 <- (y1 - y3)/(x1 - x3)
b3 <- y1 - m3*x1
hm_col_scale<-
     colorRampPalette(c("black","red"))(2)
M<-
matrix(0,ncol = 1000,nrow = 1000)
for(i in 1:1000){
     for(j in 1:1000){
           M[i,j]<- (i - m2*j - b2) > 0 & (i - m1*j - b1) > 0 & (i -m3*j - b3) <= 0
     }
}
image(M,
      col = hm_col_scale,
      zlim=c(0, 1))