SQLite Queries con el SIACON
Antonio Aguilar
10 de octubre de 2018
El Sistema de Información Agroalimentaria de Consulta (SIACON) es una base de datos publicada por el Servicio de Información Agroalimentaria y Pesquera (SIAP). La última versión fue liberada el 14 de junio de 2018. Esta base de datos contiene información sobre la producción agrícola y pecuaria de México de 1980 a 2017. La versión de descarga contiene una base de datos en SQLite (.db). Para realizar consultas en bases SQLite, R cuenta con los paquetes RSQLite y DBI.
En este documento R Markdown se realizarán consultas a la base de datos SIACON para extraer la producción nacional de frijol, así como la producción en los estados de Zacatecas, Durango y San Luis Potosí.
Establecer la conexión con SIACON
Una vez que se han cargado los paquetes RSQLite y DBI, la conexión con la base de datos se realiza como sigue:
# siacon.db se encuentra en el directorio de trabajo
con = dbConnect(SQLite(), dbname="siacon.db")Lista de tablas disponibles en siacon.db
(dbListTables(con))## [1] "agro_cierre_estatal" "agro_cierre_municipal"
## [3] "cat_ciclo" "cat_clasificacion"
## [5] "cat_cultivo" "cat_especie"
## [7] "cat_especiepesquero" "cat_estado"
## [9] "cat_generico" "cat_gruponatural"
## [11] "cat_litoral" "cat_modalidad"
## [13] "cat_municipio" "cat_producto"
## [15] "cat_tipo_produccion" "cat_tipoagricultura"
## [17] "cat_unidadmedida" "cat_variedad"
## [19] "cultivo_generico" "peco_cierre_estatal"
## [21] "peco_cierre_municipal" "pes_catalogo"
## [23] "pes_elemento" "peso_cierre_estatal"Tabla agro_cierre_estatal
Dado que se consultarán datos a nivel estatal, se usará la tabla agro_cierre_estatal. Sin embargo, para no cargar toda esa tabla, primero se abrirá la tabla cat_cultivo, para obtener el identificador que corresponde a frijol.
# Con dbGetQuery se obtiene la tabla seleccionada como data frame
cultivos <- dbGetQuery(con, "SELECT * FROM cat_cultivo") # Frijol = 6840000
head(cultivos,3)## idcultivo idunidadmedida nomcultivo cultivo
## 1 1000000 200206 Viveros de manzana Viveros de manzana
## 2 1100000 200206 Viveros de aguacate Viveros de aguacate
## 3 1200000 200206 Viveros de nuez Viveros de nuezPara encontrar el idcultivo que corresponde al frijol, se filtra el data frame de cultivos usando la columna nomcultivo.
cultivos[grep("Frijol", cultivos$nomcultivo),]## idcultivo idunidadmedida nomcultivo cultivo
## 191 6840000 200201 Frijol Frijol
## 193 6855000 200201 Frijol x pelón Frijol x pelon
## 455 15070000 200201 Frijol forrajero Frijol forrajeroEl idcultivo que corresponde al frijol es 6840000. Se seleccionan los registros bajo esa id en la tabla agro_cierre_estatal, que es la tabla que contiene los datos por estado. Esto se realiza como una consulta en SQL.
db <- dbGetQuery(con, "SELECT * FROM agro_cierre_estatal WHERE idcultivo = 6840000")
head(db)## anio idestado idciclo idmodalidad idcultivo idvariedad idunidadmedida
## 1 2016 1 2 1 6840000 6841000 200201
## 2 2016 1 2 2 6840000 6841000 200201
## 3 2016 2 2 1 6840000 6842100 200201
## 4 2016 3 1 1 6840000 6842500 200201
## 5 2016 3 2 1 6840000 6842500 200201
## 6 2016 4 1 2 6840000 6841400 200201
## idtipoagricultura idtipoproduccion sembrada cosechada siniestrada
## 1 3004 10010102 880 880.0 0.0
## 2 3004 10010102 9007 8387.0 620.0
## 3 3004 10010102 6 6.0 0.0
## 4 3004 10010102 320 320.0 0.0
## 5 3004 10010102 869 869.0 0.0
## 6 3004 10010102 2581 2564.5 16.5
## volumenproduccion valorproduccion rendimiento preciomediorural
## 1 1778.20 11227167 2.021 6313.782
## 2 2926.35 15975681 0.349 5459.252
## 3 12.00 72000 2.000 6000.000
## 4 289.19 2780367 0.904 9614.327
## 5 870.90 12907661 1.002 14821.059
## 6 1562.89 25321965 0.609 16202.013Primero, se verifica si la tabla db contiene registros a nivel nacional, bajo una idestado propia. Para ello se consulta la tabla cat_estado.
estados <- dbGetQuery(con, "SELECT * FROM cat_estado") # Estado
with(estados, unique(nomestado))## [1] "Aguascalientes" "Baja California"
## [3] "Baja California Sur" "Campeche"
## [5] "Coahuila" "Colima"
## [7] "Chiapas" "Chihuahua"
## [9] "Ciudad de México / D.F." "Durango"
## [11] "Guanajuato" "Guerrero"
## [13] "Hidalgo" "Jalisco"
## [15] "México" "Michoacán"
## [17] "Morelos" "Nayarit"
## [19] "Nuevo León" "Oaxaca"
## [21] "Puebla" "Querétaro"
## [23] "Quintana Roo" "San Luis Potosí"
## [25] "Sinaloa" "Sonora"
## [27] "Tabasco" "Tamaulipas"
## [29] "Tlaxcala" "Veracruz"
## [31] "Yucatán" "Zacatecas"Para conocer los valores de los indicatores de la tabla db, es necesario usar las tablas:
* cat_ciclo
* cat_modalidad
* cat_variedad
* cat_unidadmedida
* cat_tipoagricultura
* cat_tipo_produccion
Tabla cat_ciclo
En este punto es necesario consultar los idciclo presentes en la tabla db:
db %>%
group_by(idciclo) %>%
summarise(idciclo_count = n())## # A tibble: 2 x 2
## idciclo idciclo_count
## <int> <int>
## 1 1 2572
## 2 2 4163Para obtener los nombres que corresponden a estos idciclo, se realiza la siguiente consulta:
(ciclo <- dbGetQuery(con, "SELECT * FROM cat_ciclo WHERE idciclo IN (1,2)"))## idciclo nomciclo
## 1 1 Otoño-Invierno
## 2 2 Primavera-VeranoTabla cat_modalidad
Primero se investiga las modalidades presentes en db:
db %>%
group_by(idmodalidad) %>%
summarise(idmod_count=n())## # A tibble: 2 x 2
## idmodalidad idmod_count
## <int> <int>
## 1 1 3276
## 2 2 3459Se cuentan con dos modalidades que corresponden a:
(modalidad <- dbGetQuery(con, "SELECT * FROM cat_modalidad WHERE idmodalidad IN (1,2)"))## idmodalidad nommodalidad
## 1 1 Riego
## 2 2 TemporalTabla cat_variedad
En esta tabla se obtienen las variedades que corresponden a frijol:
variedades <- dbGetQuery(con, "SELECT * FROM cat_variedad WHERE idcultivo = 6840000")
variedades <- variedades %>%
select(-idunidadmedida,-idcultivo)
variedades$nomvariedad <- gsub(variedades$nomvariedad,pattern = "Frijol ", replacement = "")
variedades$nomvariedad <- toupper(variedades$nomvariedad)
variedades## idvariedad idgrupo nomvariedad variedad
## 1 5335100 10001060 NEGRO MICHOACÁN Frijol negro michoacan
## 2 6840100 10001060 S/CLASIFICAR Frijol s/clasificar
## 3 6840200 10001060 ALUBIA Frijol alubia
## 4 6840300 10001060 AYOCOTE Frijol ayocote
## 5 6840400 10001060 AZUFRADO Frijol azufrado
## 6 6840500 10001060 BAYO Frijol bayo
## 7 6840600 10001060 BAYO BERRENDO Frijol bayo berrendo
## 8 6840700 10001060 CACAHUATE Frijol cacahuate
## 9 6840800 10001060 CANARIO Frijol canario
## 10 6840900 10001060 FLOR DE JUNIO Frijol flor de junio
## 11 6841000 10001060 FLOR DE MAYO Frijol flor de mayo
## 12 6841100 10001060 GARBANCILLO Frijol garbancillo
## 13 6841200 10001060 MANZANO Frijol manzano
## 14 6841300 10001060 MAYOCOBA Frijol mayocoba
## 15 6841400 10001060 NEGRO JAMAPA Frijol negro jamapa
## 16 6841500 10001060 NEGRO QUERÉTARO Frijol negro queretaro
## 17 6841600 10001060 NEGRO SAN LUIS Frijol negro san luis
## 18 6841700 10001060 NEGRO VERACRUZ Frijol negro veracruz
## 19 6841800 10001060 NEGRO ZACATECAS Frijol negro zacatecas
## 20 6841900 10001060 ORGÁNICO Frijol organico
## 21 6842000 10001060 OJO DE CABRA Frijol ojo de cabra
## 22 6842100 10001060 OTROS CLAROS Frijol otros claros
## 23 6842200 10001060 OTROS DE COLOR Frijol otros de color
## 24 6842300 10001060 OTROS NEGROS Frijol otros negros
## 25 6842400 10001060 PERUANO Frijol peruano
## 26 6842500 10001060 PINTO NACIONAL Frijol pinto nacional
## 27 6842600 10001060 PINTO AMERICANO Frijol pinto americano
## 28 6842800 10001060 MARCELA Frijol marcela
## 29 6842900 10001060 PINTO SALTILLO Frijol pinto saltilloTabla cat_unidadmedida
Consultar las idunidadmedida presentes en db:
db %>%
group_by(idunidadmedida) %>%
summarise(idunidad_count=n())## # A tibble: 1 x 2
## idunidadmedida idunidad_count
## <int> <int>
## 1 200201 6735Se observa que solo hay una unidad de medida.
(unidad <- dbGetQuery(con, "SELECT * FROM cat_unidadmedida WHERE idunidadmedida IN (200201)"))## idunidadmedida nomunidad nomcorto
## 1 200201 Toneladas (ton)Tabla cat_tipoagricultura:
Se consultan los idtipoagricultura presentes en db:
db %>%
group_by(idtipoagricultura) %>%
summarise(idtipoagr_count=n())## # A tibble: 1 x 2
## idtipoagricultura idtipoagr_count
## <int> <int>
## 1 3004 6735Se observa que solo hay un tipo de agricultura en la producción de frijol.
(tipo_agr <- dbGetQuery(con, "SELECT * FROM cat_tipoagricultura WHERE idtipoagricultura IN (3004)"))## idtipoagricultura nomtipoagricultura
## 1 3004 Cielo abiertoTabla cat_tipo_produccion:
Se consultan los idtipoproduccion presentes en db:
db %>%
group_by(idtipoproduccion) %>%
summarise(idtipoprod_count=n())## # A tibble: 2 x 2
## idtipoproduccion idtipoprod_count
## <int> <int>
## 1 10010101 1
## 2 10010102 6734Se encuentran dos idtipoproduccion que corresponden a:
(tipo_prod <- dbGetQuery(con, "SELECT * FROM cat_tipo_produccion WHERE idtipoproduccion in (10010101, 10010102)"))## idtipoproduccion nomtipoproduccion
## 1 10010101 Orgánico
## 2 10010102 ConvencionalTabla agro_cierre_municipal
Información a nivel municipal se encuentra disponible desde el 2003 en la tabla agro_cierre_municipal
cultivos <- dbGetQuery(con, "SELECT * FROM agro_cierre_municipal")
head(cultivos)## anio idestado idddr idcader idmunicipio idciclo idmodalidad idcultivo
## 1 2016 1 1 1 1 1 1 5490000
## 2 2016 1 1 1 1 1 1 5900000
## 3 2016 1 1 1 1 1 1 9090000
## 4 2016 1 1 1 1 1 1 15050000
## 5 2016 1 1 1 1 2 1 5490000
## 6 2016 1 1 1 1 2 1 7330000
## idvariedad idunidadmedida idtipoagricultura idtipoproduccion sembrada
## 1 5490100 200201 3004 10010102 765
## 2 5900100 200201 3004 10010102 35
## 3 9090100 200201 3004 10010102 611
## 4 15050800 200201 3004 10010102 692
## 5 5490100 200201 3004 10010102 28
## 6 7330600 200201 3004 10010102 3
## cosechada siniestrada volumenproduccion valorproduccion rendimiento
## 1 765 0 20865.75 12424510.8 27.275
## 2 35 0 945.00 538650.0 27.000
## 3 611 0 20446.00 12127953.8 33.463
## 4 692 0 34912.70 20987769.6 50.452
## 5 28 0 746.40 407467.2 26.657
## 6 3 0 69.00 138000.0 23.000
## preciomediorural
## 1 595.45
## 2 570.00
## 3 593.17
## 4 601.15
## 5 545.91
## 6 2000.00Los años pueden consultarse como:
sort(unique(cultivos$anio))## [1] 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016
## [15] 2017Cierre de conexión
Una vez que las consultas a las tablas se terminan, se cierra la conexión a la base de datos.
dbDisconnect(con)Inspección de la tabla db
Se redefine la tabla db, para ello se deja fuera las columnas reduntantes, así como la producción orgánica, que solo tiene un registro.
db <- db %>% # Nueva definición de db
filter(idtipoproduccion == 10010102) %>% # Dejar afuera producción orgánica
select(-idcultivo, -idunidadmedida, -idtipoagricultura, -idtipoproduccion) # Dejar fuera id de cultivo, unidad de medida, tipo de agricultura y tipo de producción
db <- inner_join(db, variedades[,c("idvariedad", "nomvariedad")], by = "idvariedad")
head(db)## anio idestado idciclo idmodalidad idvariedad sembrada cosechada
## 1 2016 1 2 1 6841000 880 880.0
## 2 2016 1 2 2 6841000 9007 8387.0
## 3 2016 2 2 1 6842100 6 6.0
## 4 2016 3 1 1 6842500 320 320.0
## 5 2016 3 2 1 6842500 869 869.0
## 6 2016 4 1 2 6841400 2581 2564.5
## siniestrada volumenproduccion valorproduccion rendimiento
## 1 0.0 1778.20 11227167 2.021
## 2 620.0 2926.35 15975681 0.349
## 3 0.0 12.00 72000 2.000
## 4 0.0 289.19 2780367 0.904
## 5 0.0 870.90 12907661 1.002
## 6 16.5 1562.89 25321965 0.609
## preciomediorural nomvariedad
## 1 6313.782 FLOR DE MAYO
## 2 5459.252 FLOR DE MAYO
## 3 6000.000 OTROS CLAROS
## 4 9614.327 PINTO NACIONAL
## 5 14821.059 PINTO NACIONAL
## 6 16202.013 NEGRO JAMAPACon glimpse se puede verificar la clase de cada columna en y se aprecian los primeros registros.
glimpse(db)## Observations: 6,734
## Variables: 13
## $ anio <int> 2016, 2016, 2016, 2016, 2016, 2016, 2016, 20...
## $ idestado <int> 1, 1, 2, 3, 3, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5,...
## $ idciclo <int> 2, 2, 2, 1, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2,...
## $ idmodalidad <int> 1, 2, 1, 1, 1, 2, 2, 1, 1, 1, 2, 2, 2, 2, 2,...
## $ idvariedad <int> 6841000, 6841000, 6842100, 6842500, 6842500,...
## $ sembrada <dbl> 880.0, 9007.0, 6.0, 320.0, 869.0, 2581.0, 13...
## $ cosechada <dbl> 880.0, 8387.0, 6.0, 320.0, 869.0, 2564.5, 13...
## $ siniestrada <dbl> 0.0, 620.0, 0.0, 0.0, 0.0, 16.5, 0.0, 0.0, 0...
## $ volumenproduccion <dbl> 1778.20, 2926.35, 12.00, 289.19, 870.90, 156...
## $ valorproduccion <dbl> 11227167.2, 15975681.2, 72000.0, 2780367.4, ...
## $ rendimiento <dbl> 2.021, 0.349, 2.000, 0.904, 1.002, 0.609, 0....
## $ preciomediorural <dbl> 6313.782, 5459.252, 6000.000, 9614.327, 1482...
## $ nomvariedad <chr> "FLOR DE MAYO", "FLOR DE MAYO", "OTROS CLARO...Encontrar NA's en la tabla db:
db %>%
map_df(function(x) sum(is.na(x))) %>%
gather(feature, num_nulls)## # A tibble: 13 x 2
## feature num_nulls
## <chr> <int>
## 1 anio 0
## 2 idestado 0
## 3 idciclo 0
## 4 idmodalidad 0
## 5 idvariedad 0
## 6 sembrada 0
## 7 cosechada 0
## 8 siniestrada 0
## 9 volumenproduccion 0
## 10 valorproduccion 0
## 11 rendimiento 0
## 12 preciomediorural 0
## 13 nomvariedad 0Resumir la producción nacional y por estados
produccion <- db %>%
select(anio, idestado, volumenproduccion) %>%
group_by(anio, idestado) %>%
summarise(vol_prod=sum(volumenproduccion))
produccion <- cast(produccion, anio~idestado, value = "vol_prod")
produccion <- data.frame(produccion, row.names = "anio")
produccion$Nac <- rowSums(produccion, na.rm=T)
produccion$Reg <- rowSums(produccion[,c("X10", "X24", "X32")])
produccion %>%
select(Nac, Reg, X10, X24, X32) %>%
head()## Nac Reg X10 X24 X32
## 1980 945358 258482 141626 10021 106835
## 1981 1325034 473892 146769 29969 297154
## 1982 971697 276817 48783 26303 201731
## 1983 1279325 594983 184224 46963 363796
## 1984 923415 378212 66729 16383 295100
## 1985 907033 459357 140496 16658 302203produccion %>%
select(Nac, Reg) %>%
ggplot(data=., aes(x=as.numeric(rownames(produccion)))) +
geom_line(aes(y=Nac/1000, colour="Nac"), lwd=0.8) +
geom_line(aes(y=Reg/1000, colour="Reg"), lwd=0.8) +
scale_x_continuous(breaks = pretty(as.numeric(rownames(produccion)), n = 10)) +
ylim(0,1800) +
ylab("Producción en miles de toneladas") +
xlab("Años") +
ggtitle("Producción nacional y regional de frijol (Dur, SLP y Zac)") +
theme_bw() +
scale_color_manual(labels=c("Nac", "Reg"), values=c("Nac"="chocolate3", "Reg"="seagreen4")) +
theme(axis.text = element_text(size=9),
axis.title = element_text(size=9),
axis.line = element_line(colour = "gray"),
legend.text = element_text(size=),
legend.key = element_rect(fill = "transparent"),
legend.position = c(0.95,0.9),
legend.background = element_rect(fill = "transparent"),
legend.title = element_blank(),
panel.border = element_blank(),
panel.grid.minor = element_blank(),
panel.background = element_blank(),
plot.title = element_text(size=10, color="Black", face="bold")) 
produccion %>%
select(X10, X24, X32) %>%
ggplot(data=., aes(x=as.numeric(rownames(produccion)))) +
geom_line(aes(y=X10/1000, colour="Dur"), lwd=0.8) +
geom_line(aes(y=X24/1000, colour="SLP"), lwd=0.8) +
geom_line(aes(y=X32/1000, colour="Zac"), lwd=0.8) +
scale_x_continuous(breaks = pretty(as.numeric(rownames(produccion)), n = 10)) +
ylim(0,650) +
ylab("Producción en miles de toneladas") +
xlab("Años") +
ggtitle("Producción estatal Durango, San Luis Potosí y Zacatecas") +
theme_bw() +
scale_color_manual(labels=c("Durango","SLP","Zacatecas"), values = c("Dur"="chocolate3","SLP"="purple4","Zac"="seagreen4")) +
theme(axis.text = element_text(size=9),
axis.title = element_text(size=9),
axis.line = element_line(colour = "gray"),
legend.text = element_text(size=9),
legend.key = element_rect(fill = "transparent"),
legend.position = c(0.8,0.9),
legend.direction = "horizontal",
legend.background = element_rect(fill = "transparent"),
legend.title = element_blank(),
panel.border = element_blank(),
panel.grid.minor = element_blank(),
panel.background = element_blank(),
plot.title = element_text(size=10, color="Black", face="bold")) 
La producción se frijol negro es el estado de Zacatecas es significativa. Ahora se filtrará la producción de este estado por variedad.
db %>%
filter(idestado==32, anio > 2003) %>%
group_by(anio, nomvariedad) %>%
summarise(vol_prod = sum(volumenproduccion)) %>%
ggplot(data=.) +
geom_boxplot(aes(x=as.factor(nomvariedad), y= vol_prod/1000)) + # Gráfica en miles de toneladas
coord_flip()+
ylab("Producción en miles de toneladas") +
xlab("Variedades") +
ggtitle("Producción de frijol en Zacatecas por variedad, 2003-2017") +
theme_bw() +
theme(axis.text = element_text(size=9),
axis.title = element_text(size=9),
axis.line = element_line(colour = "gray"),
panel.border = element_blank(),
panel.grid.minor = element_blank(),
panel.background = element_blank(),
plot.title = element_text(size=10, color="Black", face="bold"))
zac_prod <- db %>%
filter(idestado==32, anio > 2003) %>%
group_by(anio, nomvariedad) %>% # Se debe agrupar, consuderando que hay dos modalidades y dos ciclos
summarise(vol_pro=sum(volumenproduccion))
zac_prod <- cast(zac_prod, anio~nomvariedad, value = "vol_pro")
colnames(zac_prod) <- gsub(x = colnames(zac_prod), pattern = " ", replacement = "_")
head(zac_prod)## anio ALUBIA BAYO CANARIO FLOR_DE_JUNIO FLOR_DE_MAYO MARCELA
## 1 2004 NA 11623.71 NA 67843.96 104155.00 NA
## 2 2005 NA 5569.12 NA 25354.74 43649.26 NA
## 3 2006 NA 19260.00 1740 76008.36 91507.22 NA
## 4 2007 NA 13952.26 1415 37706.64 48222.49 NA
## 5 2008 NA 21588.70 NA 51616.67 66450.47 NA
## 6 2009 9 18232.45 4122 29726.34 39187.49 NA
## NEGRO_SAN_LUIS NEGRO_ZACATECAS OTROS_CLAROS OTROS_NEGROS PERUANO
## 1 NA 174362.19 NA NA NA
## 2 95262.27 1849.00 NA 3164.4 NA
## 3 228302.00 5119.00 NA NA 76.0
## 4 131093.00 2772.00 NA NA 1.0
## 5 107306.00 2515.40 NA NA 20.8
## 6 107592.00 1599.14 7324 NA 7.8
## PINTO_AMERICANO PINTO_NACIONAL PINTO_SALTILLO S/CLASIFICAR
## 1 NA NA NA 6200.00
## 2 NA 27.00 NA 648.16
## 3 NA NA NA 2166.95
## 4 NA 58.95 NA 1906.35
## 5 NA NA NA 2333.60
## 6 NA 55045.00 NA 1816.75zac_prod %>%
select(anio,NEGRO_SAN_LUIS, PINTO_SALTILLO, FLOR_DE_MAYO,FLOR_DE_JUNIO) %>%
ggplot(data=., aes(x=anio)) +
geom_line(aes(y=NEGRO_SAN_LUIS/1000, colour="NSL"), lwd=0.8) +
geom_line(aes(y=PINTO_SALTILLO/1000, colour="PS"), lwd=0.8) +
geom_line(aes(y=FLOR_DE_MAYO/1000, colour="FM"), lwd=0.8) +
geom_line(aes(y=FLOR_DE_JUNIO/1000, colour="FJ"), lwd=0.8) +
scale_x_continuous(breaks = pretty(zac_prod$anio, n = 10)) +
coord_cartesian(ylim=c(0,250)) +
ylab("Producción en miles de toneladas") +
xlab("Años") +
ggtitle("Principales variedades de frijol producido en Zacatecas, 2003-2017") +
theme_bw() +
scale_color_manual(labels=c("NSL"="N San Luis","PS"="P Saltillo","FM"="F Mayo","FJ"="F Junio"),
values = c("NSL"="black","PS"="chocolate","FM"="red2", "FJ"="seagreen4")) +
theme(axis.text = element_text(size=9),
axis.title = element_text(size=9),
axis.line = element_line(colour = "gray"),
legend.text = element_text(size=9),
legend.key = element_rect(fill = "transparent"),
legend.position = c(0.75,0.9),
legend.direction = "horizontal",
legend.background = element_rect(fill = "transparent"),
legend.title = element_blank(),
panel.border = element_blank(),
panel.grid.minor = element_blank(),
panel.background = element_blank(),
plot.title = element_text(size=10, color="Black", face="bold")) ## Warning: Removed 1 rows containing missing values (geom_path).## Warning: Removed 8 rows containing missing values (geom_path).