Machine Learning Basics – Gradient Boosting & XGBoost — R-bloggers

(This article was first published on Shirin’s playgRound, and kindly contributed to R-bloggers) In a recent video, I covered Random Forests and Neural Nets as part of the codecentric.ai Bootcamp. In the most recent video, I covered Gradient Boosting and XGBoost. You can find the video on YouTube and the slides on slides.com. Both are…

via Machine Learning Basics – Gradient Boosting & XGBoost — R-bloggers

Anuncios

TRIMESTRALIZACIÓN DEL VAB INDUSTRIAL DE CANTABRIA

Lectura y tratamiento de datos

En indust_m1995.xls están los datos de los indicadores de frecuencia mensual, en indust_T1995.xls los indicadores de frecuencia trimestral, y en CRAC.xls los datos de la serie anual de la contabilidad regional de Cantabria que se quiere trimestralizar.

Datos:

indust_T1995

indust_m1995

CRAC

“`{r}
# Datos mensuales
excel1 <-read.csv(“indust_m1995.csv”, sep=”;”,header = TRUE)
excel3 <-read.csv(“indust_T1995.csv”, sep=”;”,header = TRUE)
excel4<-read.csv(“CRAC.csv”, sep=”;”,header = TRUE)
“`

Creamos objetos ts de las series de Exportaciones (EXPORT), Ocupados (OCUP) en la industria (cifras de la seguridad social), tendencia de la Cartera de Pedidos (TPEDIDOS), Indice de Precios Industriales (IPRI), Indice de producción Industrial (IPI), y Ocupados (EPA) en la Industria en la Encuesta de Población Activa.

“`{r}
EXPORT=ts(na.omit(excel1$EXPORT),frequency = 12,start=c(1995,1))
OCUP=ts(na.omit(excel1$OCUP),frequency = 12,start=c(1995,1))
TPEDIDO=ts(na.omit(excel1$TPEDIDO),frequency = 12,start=c(1995,1))
IPRI=ts(na.omit(excel1$IPRI),frequency = 12,start=c(1995,1))
IPI=ts(na.omit(excel1$IPI),frequency = 12,start=c(1995,1))
EPA=ts(na.omit(excel3$EPA),frequency = 4,start=c(1995,1))
“`

Convertimos las series de indicadores de frecuencia mensual a trimestral:

“`{r}
EXPORT=aggregate.ts(EXPORT,nfrequency = 4,mean)
OCUP=aggregate.ts(OCUP,nfrequency = 4,mean)
TPEDIDO=aggregate.ts(TPEDIDO,nfrequency = 4,mean)
IPRI=aggregate.ts(IPRI,nfrequency = 4,mean)
IPI=aggregate.ts(IPI,nfrequency = 4,mean)
“`

Deflactamos la serie de exportaciones por el índice de precios industriales.

“`{r}
EXPORT=EXPORT*100/IPRI
“`

Desestacionalizamos

Obtenemos las series de tendencia, estacionalidad con descomponer, con modelo aditivo.

“`{r}
library(descomponer)
#EXPORT
EXPORT.s=ts(descomponer(EXPORT,4,1)$datos$TDST,frequency=4,start=c(1995,1))
EXPORT.sa=ts(descomponer(EXPORT,4,1)$datos$ST,frequency=4,start=c(1995,1))
EXPORT.t=ts(descomponer(EXPORT,4,1)$datos$TD,frequency=4,start=c(1995,1))
plot(EXPORT)
lines(EXPORT.s,col=”red”)
lines(EXPORT.t,col=”blue”)
#OCUP
OCUP.s=ts(descomponer(OCUP,4,1)$datos$TDST,frequency=4,start=c(1995,1))
OCUP.sa=ts(descomponer(OCUP,4,1)$datos$ST,frequency=4,start=c(1995,1))
OCUP.t=ts(descomponer(OCUP,4,1)$datos$TD,frequency=4,start=c(1995,1))
plot(OCUP)
lines(OCUP.s,col=”red”)
lines(OCUP.t,col=”blue”)
#TPEDIDO
TPEDIDO.s=ts(descomponer(TPEDIDO,4,1)$datos$TDST,frequency=4,start=c(1995,1))
TPEDIDO.sa=ts(descomponer(TPEDIDO,4,1)$datos$ST,frequency=4,start=c(1995,1))
TPEDIDO.t=ts(descomponer(TPEDIDO,4,1)$datos$TD,frequency=4,start=c(1995,1))
plot(TPEDIDO)
lines(TPEDIDO.s,col=”red”)
lines(TPEDIDO.t,col=”blue”)
#IPI
IPI.s=ts(descomponer(IPI,4,1)$datos$TDST,frequency=4,start=c(1995,1))
IPI.sa=ts(descomponer(IPI,4,1)$datos$ST,frequency=4,start=c(1995,1))
IPI.t=ts(descomponer(IPI,4,1)$datos$TD,frequency=4,start=c(1995,1))
plot(IPI)
lines(IPI.s,col=”red”)
lines(IPI.t,col=”blue”)
#EPA
EPA.s=ts(descomponer(EPA,4,1)$datos$TDST,frequency=4,start=c(1995,1))
EPA.sa=ts(descomponer(EPA,4,1)$datos$ST,frequency=4,start=c(1995,1))
EPA.t=ts(descomponer(EPA,4,1)$datos$TD,frequency=4,start=c(1995,1))
plot(EPA)
lines(EPA.s,col=”red”)
lines(EPA.t,col=”blue”)
“`

Pronósticos de la serie desestacionalizada

Con tslm de forecast pronosticamos: tendencias y estacionalidades. Prolongamos la serie original (o) con el pronostico de tendencia y estacionalidad, la serie de tendencia y estacionalidad (f) con el pronostico de tendencia y estacionalidad, y la serie de tendencia (t) con el pronóstico de tendencia.

“`{r}
library(forecast)
#EXPORT
fit1 <- auto.arima(EXPORT.sa)
fit2 <- auto.arima(EXPORT.t)
plot(forecast(fit1, h=20))
plot(forecast(fit2, h=20))
fit3 = forecast(fit1, h=20)
fit4 = forecast(fit2, h=20)
EXPORT.f=ts(c(EXPORT.s,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
EXPORT.o=ts(c(EXPORT,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
EXPORT.t=ts(c(EXPORT.t,fit4$mean),frequency = 4,start=c(1995,1))
plot(EXPORT.o)
lines(EXPORT.f,col=”blue”)
lines(EXPORT.t,col=”red”)
#OCUP
fit1 <- auto.arima(OCUP.sa)
fit2 <- auto.arima(OCUP.t)
plot(forecast(fit1, h=20))
plot(forecast(fit2, h=20))
fit3 = forecast(fit1, h=20)
fit4 = forecast(fit2, h=20)
OCUP.f=ts(c(OCUP.s,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
OCUP.o=ts(c(OCUP,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
OCUP.t=ts(c(OCUP.t,fit4$mean),frequency = 4,start=c(1995,1))
plot(OCUP.o)
lines(OCUP.f,col=”blue”)
lines(OCUP.t,col=”red”)
# TPEDIDO
fit1 <- auto.arima(TPEDIDO.sa)
fit2 <- auto.arima(TPEDIDO.t)
plot(forecast(fit1, h=20))
fit3 = forecast(fit1, h=20)
fit4 = forecast(fit2, h=20)
TPEDIDO.f=ts(c(TPEDIDO.s,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
TPEDIDO.o=ts(c(TPEDIDO,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
TPEDIDO.t=ts(c(TPEDIDO.t,fit4$mean),frequency = 4,start=c(1995,1))
plot(TPEDIDO)
lines(TPEDIDO.f,col=”blue”)
lines(TPEDIDO.t,col=”red”)
#IPI
fit1 <- auto.arima(IPI.sa)
fit2 <- auto.arima(IPI.t)
plot(forecast(fit1, h=20))
fit3 = forecast(fit1, h=20)
fit4 = forecast(fit2, h=20)
IPI.f=ts(c(IPI.s,c(fit3$mean)+c(fit4$mean)),frequency = 4,start=c(1995,1))
IPI.o=ts(c(IPI,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
IPI.t=ts(c(IPI.t,fit4$mean),frequency = 4,start=c(1995,1))
plot(IPI)
lines(IPI.f,col=”blue”)
lines(IPI.t,col=”red”)
# EPA
fit1 <- auto.arima(EPA.sa)
fit2 <- auto.arima(EPA.t)
plot(forecast(fit1, h=20))
fit3 = forecast(fit1, h=20)
fit4 = forecast(fit2, h=20)
EPA.f=ts(c(EPA.s,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
EPA.o=ts(c(EPA,fit3$mean+fit4$mean),frequency = 4,start=c(1995,1))
EPA.t=ts(c(EPA.t,fit4$mean),frequency = 4,start=c(1995,1))
plot(EPA)
lines(EPA.f,col=”blue”)
lines(EPA.t,col=”red”)

“`

Agregación de series

Obtenemos los indicadores en frecuencias anuales desde las frecuencias trimestrales.

“`{r}
EXPORT.A=aggregate.ts(EXPORT.f,nfrequency = 1,mean)
OCUP.A=aggregate.ts(OCUP.f,nfrequency = 1,mean)
TPEDIDO.A=aggregate.ts(TPEDIDO.f,nfrequency = 1,mean)
IPI.A=aggregate.ts(IPI.f,nfrequency = 1,mean)
EPA.A=aggregate.ts(EPA.f,nfrequency = 1,mean)
“`

Creamos varios data frame con las series anuales para el periodo 2000-2019. En cada data frame ha

“`{r}
fin=2019
# SELECCION BASE
INDUS_CAN.A=data.frame(EXPORT=window(EXPORT.A,start=c(2000),end=c(fin)),OCUP=window(OCUP.A,start=c(2000),end=c(fin)),IPI=window(IPI.A,start=c(2000),end=c(fin)))
INDUS_CAN.A.ts=ts(INDUS_CAN.A,start=c(2000))
INDUS_CAN.A.tasa=diff(INDUS_CAN.A.ts)/INDUS_CAN.A.ts
colnames(INDUS_CAN.A.tasa)=c(“EXPORT”,”OCUP”,”IPI”)
# SELECCION 1
INDUS_CAN.A.1=data.frame(EXPORT=window(EXPORT.A,start=c(2000),end=c(fin)),OCUP=window(OCUP.A,start=c(2000),end=c(fin)),TPEDIDO=window(TPEDIDO.A,start=c(2000),end=c(fin)),IPI=window(IPI.A,start=c(2000),end=c(fin)))
INDUS_CAN.A.1.ts=ts(INDUS_CAN.A.1,start=c(2000))
INDUS_CAN.A.1.tasa=diff(INDUS_CAN.A.1.ts)/INDUS_CAN.A.1.ts
colnames(INDUS_CAN.A.1.tasa)=c(“EXPORT”,”OCUP”,”TPEDIDO”,”IPI”)
# SELECCION 2
INDUS_CAN.A.2=data.frame(EXPORT=window(EXPORT.A,start=c(2000),end=c(fin)),IPI=window(IPI.A,start=c(2000),end=c(fin)),EPA=window(EPA.A,start=c(2000),end=c(fin)))
INDUS_CAN.A.2.ts=ts(INDUS_CAN.A.2,start=c(2000))
INDUS_CAN.A.2.tasa=diff(INDUS_CAN.A.2.ts)/INDUS_CAN.A.2.ts
colnames(INDUS_CAN.A.2.tasa)=c(“EXPORT”,”IPI”,”EPA”)
# SELECCION 3
INDUS_CAN.A.3=data.frame(OCUP=window(OCUP.A,start=c(2000),end=c(fin)),IPI=window(IPI.A,start=c(2000),end=c(fin)))
INDUS_CAN.A.3.ts=ts(INDUS_CAN.A.3,start=c(2000))
INDUS_CAN.A.3.tasa=diff(INDUS_CAN.A.3.ts)/INDUS_CAN.A.3.ts
colnames(INDUS_CAN.A.3.tasa)=c(“OCUP”,”IPI”)
“`

Leemos la serie anual del Indice de Volumen del VAB de la industria de la CRE y obtenemos la tasas de crecimento anual:

“`{r}
VAB=ts(na.omit(excel4$VAB2_Co),start=c(2000))
VAB.tasa=diff(VAB)/VAB
“`

Índices con Análisis factorial

Elaboramos un indicador de síntesis con la selección básica de indicadores, utilizando la librería “rela”, se correlaciona el indicador obtenido con la tasa anual del crecimiento del VAB en la industria, y se proyecta la serie anual del VAB hasta 2019.

“`{r}
library(rela)
res <- paf(as.matrix(INDUS_CAN.A.tasa))
summary(res) # Obtenemos el KMO con MSA, para determinar el número de factores
# calculamos el chi-square del Bartlett’s sphericity test, la comunalidad y
# pesos de los factores.
barplot(res$Eigenvalues[,1]) # Representamos los eigenvalues.
resv <- varimax(res$Factor.Loadings) # Rotación Varimax
print(resv)
#barplot(sort(colSums(loadings(resv)^2),decreasing=TRUE)) # Gráfico utilizando los pesos rotados.
scores <- scale(as.matrix(INDUS_CAN.A.tasa)) %*% as.matrix(resv) # Pesos factoriales, para obtener los indices tambien se puede
# Intervenciones
d13=c(rep(0,12),1,rep(0,(17-13)))
d14=c(rep(0,13),1,rep(0,(17-14)))
# Modelo 1: Estimacion MCO
indice1=ts(scores[,1],start=c(2001))
indice1.1=window(indice1,end=c(2017))
mod1=lm(VAB.tasa~0+indice1.1)
summary(mod1)
plot(mod1)
# Grafico
plot(VAB.tasa)
lines(ts(mod1$fitted,start=c(2001)),col=”red”)
# Pronósticos
indice1.1=window(indice1,start=c(2018),end=c(2019))
nuevos=data.frame(indice1.1)
VAB.tasa.fa=ts(predict(mod1,nuevos,interval=”prediction”),start=c(2018),end=c(2019))
VAB.tasa.fa
VAB.2018=window(VAB,start=c(2017),end=c(2017))*c((1+window(VAB.tasa.fa[,”fit”],start=c(2018),end=c(2018))))
VAB.2019=VAB.2018*c((1+window(VAB.tasa.fa[,”fit”],start=c(2019),end=c(2019))))
VAB1=ts(c(VAB,VAB.2018,VAB.2019),start=c(2000))
plot(VAB1)
“`

Índices con el método de Granger y Newbold

Elaboramos un indicador de síntesis para la industria utilizando el método de Granger y Newbold y la selección básica de indicadores, se correlaciona el indicador obtenido con la tasa anual del crecimiento del VAB en la industria, y se proyecta la serie anual del VAB hasta 2019 :

“`{r}

# Indice
sd1=1/sd(lm(VAB.tasa~window(INDUS_CAN.A.tasa[,”EXPORT”],end=c(2017)))$residuals)
sd2=1/sd(lm(VAB.tasa~window(INDUS_CAN.A.tasa[,”OCUP”],end=c(2017)))$residuals)
sd3=1/sd(lm(VAB.tasa~window(INDUS_CAN.A.tasa[,”IPI”],end=c(2017)))$residuals)
sdT=sd1+sd2+sd3
indice2=(sd1/sdT)*INDUS_CAN.A.tasa[,”EXPORT”]+(sd2/sdT)*INDUS_CAN.A.tasa[,”OCUP”]+(sd3/sdT)*INDUS_CAN.A.tasa[,”IPI”]
indice2=ts(indice2,start=c(2001))
indice2.1=window(indice2,end=c(2017))
# Estimacion MCO
mod2=lm(VAB.tasa~0+indice2.1)
#mod2=lm(VAB.tasa~0+indice2.1)
summary(mod2)
plot(mod2)
# Grafico
plot(VAB.tasa)
lines(ts(mod2$fitted,start=c(2001)),col=”red”)
# Pronósticos
indice2.1=window(indice2,start=c(2018),end=c(2019))
nuevos=data.frame(indice2.1)
VAB.tasa.f=ts(predict(mod2,nuevos,interval=”prediction”),start=c(2018),end=c(2019))
VAB.tasa.f
VAB.2018=window(VAB,start=c(2017),end=c(2017))*c((1+window(VAB.tasa.f[,”fit”],start=c(2018),end=c(2018))))
VAB.2019=VAB.2018*c((1+window(VAB.tasa.f[,”fit”],start=c(2019),end=c(2019))))
VAB2=ts(c(VAB,VAB.2018,VAB.2019),start=c(2000))
plot(VAB2)
“`

Trimestralizar

Trimestralizamos la serie original con el método de Granger y Newbold, utilizando la librería tempdisagg, y el método de Chow-Lee.

“`{r}

library(tempdisagg)
indice.o=(sd1/sdT)*scale(window(EXPORT.o,start=c(2000,1),end=c(2018,1)))+(sd2/sdT)*scale(window(OCUP.o,start=c(2000,1),end=c(2018,1)))+(sd3/sdT)*scale(window(IPI.o,start=c(2000,1),end=c(2018,1)))
indice.o=ts(indice.o,frequency=4,start=c(2000,1))
mod.o=td(VAB1~indice.o,method = “chow-lin-maxlog”)
summary(mod.o)
plot(mod.o)
resultado.o=predict(mod.o)
resultado.o
plot(resultado.o)
“`

Trimestralizamos la serie de tendencia, y la serie de tendencia y estacionalidad, utilizando la libreria tempdisagg, y el método de Chow-Lee, representamos los resultados finales:

“`{r}
indice.f=(sd1/sdT)*scale(window(EXPORT.f,start=c(2000,1),end=c(2018,1)))+(sd2/sdT)*scale(window(OCUP.f,start=c(2000,1),end=c(2018,1)))+(sd3/sdT)*scale(window(IPI.f,start=c(2000,1),end=c(2018,1)))
indice.f=ts(indice.f,frequency=4,start=c(2000,1))
mod.f=td(VAB1~indice.f,method = “chow-lin-maxlog”)
summary(mod.f)
plot(mod.f)
indice.T=(1/sd1)/(1/sdT)*scale(window(EXPORT.t,start=c(2000,1),end=c(2018,1)))+(1/sd2)/(1/sdT)*scale(window(OCUP.t,start=c(2000,1),end=c(2018,1)))+(1/sd3)/(1/sdT)*scale(window(IPI.t,start=c(2000,1),end=c(2018,1)))
indice.T=ts(indice.T,frequency=4,start=c(2000,1))
mod.T=td(VAB1~indice.T,method = “chow-lin-maxlog”)
summary(mod.T)
plot(mod.T)
resultado.f=predict(mod.f)
resultado.f
resultado.T=predict(mod.T)
resultado.T
plot(resultado.o)
lines(resultado.T,col=”red”)
lines(resultado.f,col=”blue”)
“`

Resutados en word

Trimestralizacion industria

En R pub:

https://rpubs.com/PacoParra/443765

 

 

 

Machine learning logistic regression for credit modelling in R — R-bloggers

(This article was first published on R Programming – DataScience+, and kindly contributed to R-bloggers) Categories Regression Models Tags Logistic Regression R Programming ROC Machine learning logistic regressions is a widely popular method to model credit modeling. There are excellent and efficient packages in R, that can perform these types of analysis. Typically you will…

via Machine learning logistic regression for credit modelling in R — R-bloggers

5 Alternatives to the Default R Outputs for GLMs and Linear Models — R-bloggers

(This article was first published on R – Displayr, and kindly contributed to R-bloggers) The standard summary outputs from the glm and lm summary methods are a case in point. If you have been using R for as long as I have (19 or 20 years…) you will no doubt have a certain affection for…

via 5 Alternatives to the Default R Outputs for GLMs and Linear Models — R-bloggers

Recreating Michael Betancourt’s Bayesian modeling course from his online materials — R-bloggers

 

(This article was first published on Shravan Vasishth’s Slog (Statistics blog), and kindly contributed to R-bloggers)

Several people wanted to have the slides from Betancourt’s lectures at SMLP2018. It is possible to recreate most of the course from his writings:

1. Intro to probability:https://betanalpha.github.io/assets/case_studies/probability_theory.html

2. workflow:https://betanalpha.github.io/assets/case_studies/principled_bayesian_workflow.html

3. Diagnosis:https://betanalpha.github.io/assets/case_studies/divergences_and_bias.html

4. HMC: https://www.youtube.com/watch?v=jUSZboSq1zg

5. Validating inference: https://arxiv.org/abs/1804.06788 6. Calibrating inference:…

via Recreating Michael Betancourt’s Bayesian modeling course from his online materials — R-bloggers

Análisis espacial con R: Usa R como un Sistema de Información Geográfica.

El libro “Análisis espacial con R: Usa R como un Sistema de Información Geográfica” ha sido escrito por JEAN-FRANCOIS MAS, data de julio de 2018 y está bajo licencia libre CC.

Análisis espacial con R

R es una plataforma de análisis estadístico con herramientas gráficas muy avanzadas, es un referente en el análisis estadístico desde hace muchos años.

A primera vista, R puede parecer poco amigable a usuarios acostumbrados a manejar programas computacionales con menús y opciones seleccionadas con el ratón debido a que se basa en líneas de comando. Sin embargo, después de haber (fácilmente) superado este obstáculo, estos usuarios verán que el uso de pequeños guiones ”scripts” que permiten ejecutar una secuencia de operaciones, es mucho más eficiente que una larga secuencia de ”clics”, sin olvidar la reducción del riesgo de tendinitis. Permite repetir fácilmente el mismo procedimiento con datos diferentes o realizar modificaciones a una cadena de procesamiento ya implementada. Adicionalmente, reduce enormemente la posibilidad de cometer errores en una cadena de operaciones rutinarias y permite documentar el procesamiento realizado.

El presente libro se dirige a usuarios con conocimiento básico de Sistemas de Información Geográfica (SIG) que desean iniciarse en el manejo y análisis de datos espaciales en R. No requiere por lo tanto de ningún conocimiento previo de este programa pero si un conocimiento básico de los SIG. El libro pretende permitir al lector dar los primeros pasos en el manejo de R para el análisis espacial sin demasiados tropiezos. Para seguir con aplicaciones más avanzadas, existe un gran número de fuentes de información (ver anexos).

El libro “Análisis espacial con R” se organizó de la siguiente manera: en el primer capítulo se explica como instalar R y RStudio y se presentan los principales elementos de la interface RStudio. Se recomienda realizar la instalación de ambos programas para poder experimentar los códigos de los capítulos siguientes. En el segundo capítulo, se hace una iniciación al manejo básico de R. El lector con conocimiento previo de R puede pasar directamente al siguiente capítulo. En el tercer capítulo, se presenta como están estructurados los datos espaciales en R en los paquetes sf y raster, los dos paquetes que vamos utilizar a lo largo de este libro. La estructura de los datos en el paquete sp se encuentra en anexos. Este capítulo puede parecer un poco árido. De hecho se puede manejar información espacial sin entrar en los detalles de la organización de la información. Sin embargo, es importante, y ayuda mucho, conocer esta información. En el capítulo 4, se presentan algunas formas para intercambiar datos geográficos entre R y otros sistemas de manejo de información geográfica a través de procedimientos de importación / exportación entre R y datos en formato vectorial o de imagen, así como algunos métodos para convertir información entre vector y raster. En los capítulos 5 y 6, se presentan operaciones básicas de SIG, respectivamente con datos en formato vector y raster. En el capítulo 7, se muestran algunos de los numerosos análisis de tipo geoestadístico que se puede llevar a cabo con paquetes de R. En el octavo capítulo, se aborda el análisis de imágenes de satélite. En capítulo 9 muestra algunas formas de elaborar cartografía. Finalmente, el décimo capítulo introduce al lector las técnicas para hacer interactuar R con el programa SIG de código abierto Q-GIS y la plataforma de modelación espacial Dinámica EGO.

Descarga del libro

Mas, J-F., 2018, Análisis espacial con R: Usa R como un Sistema de Información Geográfica, European Scientific Institute, 114 p (+ anexos).

Enlace de descarga: eujournal.org/files/journals/1/books/JeanFrancoisMas.pdf

ISBN: 978-608-4642-66-4.