Analisis de datos historico de linea de tiempo 1929-2022

knitr::opts_chunk$set(
  echo = FALSE, message = F,
  warning = F, ft.keepnext = F, tab.topcaption = TRUE,
  ft.align = "left", tab.cap.pre = "Tabla ", tab.cap.sep = ": "
)
library(tidyverse)
library(emmeans)
library(agricolae)
library(janitor)
library(data.table)
library(cluster)
library(forcats)
library(gower)
library(officer)
library(ggpubr)
library(flextable)
library(knitr)
library(dumbbell)
library(ggplot2)
library(sf)
library(rnaturalearth)
library(rnaturalearthdata)
library(ComplexUpset)
### librarias para modelos
library(gmodels)
library(nnet)
library(emmeans)
source("D:/omar/METRIKA-GROUP/Github/agrograf/count_utils.R")
source("D:/omar/METRIKA-GROUP/Github/agrograf/upset_tools.R")
source("D:/omar/METRIKA-GROUP/Github/agrograf/upset_plot.R")
source("D:/omar/METRIKA-GROUP/Github/agrograf/bar_plot.R")
source("R/utils.R")
source("R/exsity_analysis_genebank.R")
set.seed(123)
autonum <- run_autonum(seq_id = "tab", bkm = "TC1", bkm_all = TRUE)
options(knitr.duplicate.label = "allow")
datos_time <- readxl::read_xlsx(path = "D:/omar/METRIKA-GROUP/Github/analisis_datos_cip/TimeLine/Data_Curada/curado_inventario_general/curado_inventario_1927-202_v04.xlsx") %>% as.data.frame()

## 

1 Descripción del estudio

El Estudio de Linea de Tiempo Historico contiene un total de 12000. Los datos del estudio estan comprendidos entre los años 1927 y 2022. El estudio fue llevado a cabo en Paucartambo, en la región Cusco, en los distritos de Paucartambo, Colquepata, Huancarani, Challabamba, Caicay y Kosñipata.

Debajo una tabla descriptiva con las observaciones estratificadas por año:

res <- datos_time %>% tabyl(ADM3_Name, cu_date_collection)
adm3 <- c(res$ADM3_Name,"Total")
res <- res %>% adorn_totals("row") %>% select(-ADM3_Name) %>% mutate(Total = rowSums(.))
names(res) <- paste0("y", names(res))
res <- res %>% add_column(ADM3_Name=adm3) %>% relocate(ADM3_Name, .before = "y1945")
names(res) <- str_remove_all(names(res),"y")
flextable(res)

ADM3_Name

1945

1947

1951

1953

1963

1970

1971

1972

1975

1982

1986

1988

1990

1991

1993

1994

1997

2001

2005

2008

2009

2013

2016

2019

2022

Total

Challabamba

0

87

0

0

0

0

0

3

12

0

42

0

0

0

0

0

0

0

0

1

2,556

0

0

491

2,935

6,127

Colquepata

0

128

0

0

0

2

0

0

29

0

272

73

0

50

4

110

0

460

0

3

0

26

8

531

2,143

3,839

Paucartambo

50

167

1

1

1

9

70

10

172

1

515

0

21

0

15

109

31

460

1

1

0

8

8

802

2,838

5,291

Total

50

382

1

1

1

11

70

13

213

1

829

73

21

50

19

219

31

920

1

5

2,556

34

16

1,824

7,916

15,257

Numero de observaciones estratificada por año y distrito

ADM3_Name

1941-1960

1961-1980

1981-2000

2001-2020

2021-actualidad

Total

Challabamba

87

15

42

3,048

2,935

6,127

Colquepata

128

31

509

1,028

2,143

3,839

Paucartambo

219

262

692

1,280

2,838

5,291

Total

434

308

1,243

5,356

7,916

15,257

Nota

Se observa entonces que los datos con mayor presencia de observaciones son Paucartambo, Challabamba y Colquepata.

Se observa que el numero de observaciones es completo a partir de 1941-1960, puesto que no existe datos de Challabamba para la decada de 1930-1940.

Teniendo en cuenta las anteriores consideraciones, filtramos los datos con esos criterios

2 Estadística descriptiva

2.1 Frecuencia de variedades unicas por periodo de años

f_period <- count_distinct(datos_time, "period_date_strat", "cu_variety_name",percentage = TRUE)
flextable(f_period)

period_date_strat

n_cu_variety_name

prcnt_cu_variety_name

1941-1960

226

8.9

1961-1980

212

8.4

1981-2000

353

13.9

2001-2020

699

27.6

2021-actualidad

1,047

41.3

2.2 Frecuencia de variedades unicas por cultivar group

f_cgroup <- count_distinct(datos_time, "cultivar_group", "cu_variety_name",percentage = TRUE)
flextable(f_cgroup)

cultivar_group

n_cu_variety_name

prcnt_cu_variety_name

Bitter landrace

138

9.0

Floury landrace

1,368

89.4

Modern variety

24

1.6

2.3 Frecuencia de variedades únicas por ADM3

f_adm3 <- count_distinct(datos_time, "ADM3_Name", "cu_variety_name",percentage = TRUE)
flextable(f_adm3 )

ADM3_Name

n_cu_variety_name

prcnt_cu_variety_name

Challabamba

645

26.6

Colquepata

812

33.4

Paucartambo

972

40.0

2.4 Frecuencia de variedades únicas por cultivar group y periodo de años

f_adm3_period <- count_distinct(datos_time, c("ADM3_Name", "period_date_strat") ,"cu_variety_name") %>% group_by(ADM3_Name) %>% mutate(prcnt_cu_variety_name= pct_fun(n_cu_variety_name))
flextable(f_adm3_period )

ADM3_Name

period_date_strat

n_cu_variety_name

prcnt_cu_variety_name

Challabamba

1941-1960

65

7.1

Challabamba

1961-1980

15

1.6

Challabamba

1981-2000

30

3.3

Challabamba

2001-2020

302

32.9

Challabamba

2021-actualidad

506

55.1

Colquepata

1941-1960

93

7.5

Colquepata

1961-1980

27

2.2

Colquepata

1981-2000

200

16.0

Colquepata

2001-2020

377

30.2

Colquepata

2021-actualidad

551

44.2

Paucartambo

1941-1960

151

9.2

Paucartambo

1961-1980

189

11.5

Paucartambo

1981-2000

260

15.8

Paucartambo

2001-2020

530

32.2

Paucartambo

2021-actualidad

517

31.4

2.5 Frecuencia de variedades unicas por adm3, cultivar group y periodo de años

f_adm3_vari_period <- count_distinct(datos_time, c("ADM3_Name", "period_date_strat","cultivar_group") ,"cu_variety_name") %>% group_by(ADM3_Name,period_date_strat) %>% mutate(pct= pct_fun(n_cu_variety_name)) %>% ungroup()

2.6 Presencia de variedades (únicas) más comunes que estan presentes en diferentes periodos

upstbl <- freq_upset(datos_time, "cu_variety_name", "period_date_strat")
upspam <- pam_upset(upstbl)
graf_ups_period_variety <- upset_plot(upspam, c("1941-1960","1961-1980","1981-2000", "2001-2020", "2021-actualidad"), c("black","green","blue","red","brown"))
graf_ups_period_variety

2.7 Listado de variedades (conteo único) más comunes en todos los periodos de años

tbl_vunica_period <- upspam %>% filter(`2021-actualidad`==1, `2001-2020`==1, `1981-2000`==1,`1961-1980`==1,`1941-1960`==1) %>% select(cu_variety_name)
vct_vunica_period <- tbl_vunica_period$cu_variety_name %>% as.character()

res <- datos_time %>%  filter(cu_variety_name %in%  vct_vunica_period)
out <- janitor::tabyl(res, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period <- out[nrow(out),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period) <- c("category","n")
top_comunes_period <- top_comunes_period %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n)) 
DT::datatable(top_comunes_period)
hbar_plot(top_comunes_period %>% slice(1:10), "category","pct",title = "Top 10 variedades más comunes en 1930-2022")

2.8 Listado de variedades unicas en el periodo 1941-1960

Esta lista contiene un conjunto de variedades que solamente están presentes en el periodo 1941-1960.

tbl_vunica_period_4160 <- upspam %>% filter(`2021-actualidad`==0, `2001-2020`==0, `1981-2000`==0,`1961-1980`==0,`1941-1960`==1) %>% select(cu_variety_name)
vunica_period_4160 <- tbl_vunica_period_4160$cu_variety_name %>% as.character()

res_4160 <- datos_time %>%  filter(cu_variety_name %in%  vunica_period_4160)

El total de variedades unicas durantel periodo 29.

out4160 <- janitor::tabyl(res_4160, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period4160 <- out4160[nrow(out4160),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period4160) <- c("category","n")
top_comunes_period_4160 <- top_comunes_period4160 %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n))

Tabla de las variedades y frecuencias de las variedades más comunes de 1941-1960

DT::datatable(top_comunes_period_4160)

2.9 Listado de variedades unicas en el periodo 1961-1980

Esta lista contiene un conjunto de variedades que solamente están presentes en el periodo 1961-1980.

tbl_vunica_period_6180 <- upspam %>% filter(`2021-actualidad`==0, `2001-2020`==0, `1981-2000`==0,`1961-1980`==1,`1941-1960`==0) %>% select(cu_variety_name)
vunica_period_6180 <- tbl_vunica_period_6180$cu_variety_name %>% as.character()

res_6180 <- datos_time %>%  filter(cu_variety_name %in%  vunica_period_6180)

El total de variedades unicas durantel periodo 42.

out6180 <- janitor::tabyl(res_6180, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period6180 <- out4160[nrow(out6180),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period6180) <- c("category","n")
top_comunes_period_6180 <- top_comunes_period6180 %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n))
DT::datatable(top_comunes_period_6180)

2.10 Listado de variedades unicas en el periodo 1981-200

Esta lista contiene un conjunto de variedades que solamente están presentes en el periodo 1981-2000.

tbl_vunica_period_8120 <- upspam %>% filter(`2021-actualidad`==0, `2001-2020`==0, `1981-2000`==1,`1961-1980`==0,`1941-1960`==0) %>% select(cu_variety_name)
vunica_period_8120 <- tbl_vunica_period_8120$cu_variety_name %>% as.character()

res_8120 <- datos_time %>%  filter(cu_variety_name %in%  vunica_period_8120)

El total de variedades unicas durantel periodo 121.

out8120 <- janitor::tabyl(res_8120, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period8120 <- out8120[nrow(out8120),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period8120) <- c("category","n")
top_comunes_period_8120 <- top_comunes_period8120 %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n))
DT::datatable(top_comunes_period_6180)

2.11 Listado de variedades unicas en el periodo 2000-2020

Esta lista contiene un conjunto de variedades que solamente están presentes en el periodo 2000-2020.

tbl_vunica_period_200020 <- upspam %>% filter(`2021-actualidad`==0, `2001-2020`==1, `1981-2000`==0,`1961-1980`==0,`1941-1960`==0) %>% select(cu_variety_name)
vunica_period_200020 <- tbl_vunica_period_200020$cu_variety_name %>% as.character()

res_200020 <- datos_time %>%  filter(cu_variety_name %in%  vunica_period_200020)

El total de variedades unicas durantel periodo 232.

out200020 <- janitor::tabyl(res_200020, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period200020 <- out200020[nrow(out200020),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period200020) <- c("category","n")
top_comunes_period_200020 <- top_comunes_period200020 %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n))
DT::datatable(top_comunes_period_200020)

2.12 Listado de variedades unicas en el periodo 2020-actualidad

Esta lista contiene un conjunto de variedades que solamente están presentes en el periodo 2020-actualidad.

tbl_vunica_period_21actual <- upspam %>% filter(`2021-actualidad`==1, `2001-2020`==0, `1981-2000`==0,`1961-1980`==0,`1941-1960`==0) %>% select(cu_variety_name)
vunica_period_21actual <- tbl_vunica_period_21actual$cu_variety_name %>% as.character()

res_21actual <- datos_time %>%  filter(cu_variety_name %in%  vunica_period_21actual)

El total de variedades unicas durantel periodo 527.

out21actual <- janitor::tabyl(res_21actual, period_date_strat ,cu_variety_name) %>% 
        janitor::adorn_totals("row")
top_comunes_period21actual <- out21actual[nrow(out21actual),] %>% t() %>% as.data.frame() %>% rownames_to_column() %>%  slice(-1)
colnames(top_comunes_period21actual) <- c("category","n")
top_comunes_period_21actual <- top_comunes_period21actual %>% mutate(n=as.numeric(n)) %>%   mutate(pct=pct_fun(n)) %>% arrange(desc(n))
DT::datatable(top_comunes_period_21actual)

2.13 Presencia de variedades (únicas) más comunes que estan presentes en diferentes distritos (Challambamba, Colquepata, Paucartambo)

upstbl <- freq_upset(datos_time, "cu_variety_name", "ADM3_Name")
upspam <- pam_upset(upstbl)
graf_pres_adm3 <- upset_plot(upspam, c("Challabamba", "Colquepata", "Paucartambo"), c("green","blue","red"))
graf_pres_adm3

hbar_freq_plot(top_comunes_adm3 %>% slice(1:10),fill_color = "lightgrey", "category","n",upper_limit = 500 ,title = "Top 10 variedades más comunes en Paucartambo-Challabamba-Colquepata")

2.14 Presencia de variedades más comunes (conteo unico) en los periodos y los distritos

datos_time <- datos_time %>% mutate(perdis=paste(ADM3_Name,period_date_strat,sep="_"))
upstbl <- freq_upset(datos_time, "cu_variety_name", "perdis")
upspam <- pam_upset(upstbl)
graf_ups_period_variety <- upset_plot(upspam, c("Challabamba_1941-1960",
"Challabamba_1961-1980","Challabamba_1981-2000","Challabamba_2001-2020",
"Challabamba_2021-actualidad","Colquepata_1941-1960",
"Colquepata_1961-1980","Colquepata_1981-2000","Colquepata_2001-2020",
"Colquepata_2021-actualidad","Paucartambo_1941-1960",
"Paucartambo_1961-1980","Paucartambo_1981-2000","Paucartambo_2001-2020",
"Paucartambo_2021-actualidad"), colors = rep("black",15))
top10_comunes_dis_period<- upspam %>% 
    adorn_totals(where="col") %>% 
    arrange(desc(Total)) %>% 
    relocate(Total, .before="Challabamba_1941-1960")  %>%  
    slice(1:10)

#top10_comunes_dis_period

flextable(top10_comunes_dis_period)

cu_variety_name

Total

Challabamba_1941-1960

Challabamba_1961-1980

Challabamba_1981-2000

Challabamba_2001-2020

Challabamba_2021-actualidad

Colquepata_1941-1960

Colquepata_1961-1980

Colquepata_1981-2000

Colquepata_2001-2020

Colquepata_2021-actualidad

Paucartambo_1941-1960

Paucartambo_1961-1980

Paucartambo_1981-2000

Paucartambo_2001-2020

Paucartambo_2021-actualidad

Huaman Uma

14

1

1

1

1

1

0

1

1

1

1

1

1

1

1

1

Puka Bole

14

1

1

1

1

1

1

0

1

1

1

1

1

1

1

1

Yana Cusi

14

1

1

1

1

1

1

0

1

1

1

1

1

1

1

1

Chiqchi Puru

13

1

1

0

1

1

1

0

1

1

1

1

1

1

1

1

Muru Qusi

13

1

1

0

1

1

1

0

1

1

1

1

1

1

1

1

Phaspa Sunchu

13

1

1

0

1

1

1

1

1

1

1

0

1

1

1

1

Puka Mama

13

1

1

1

1

1

0

0

1

1

1

1

1

1

1

1

Qequrani

13

1

0

1

1

1

1

0

1

1

1

1

1

1

1

1

Yana Bole

13

1

0

0

1

1

1

1

1

1

1

1

1

1

1

1

Muru Bole

12

1

0

0

1

1

1

0

1

1

1

1

1

1

1

1

hbar_freq_plot(top10_comunes_dis_period , "cu_variety_name","Total",title = "Top 10 variedades más comunes en 1930-2022")

2.15 Presencia de variedades que estan en conservación insitu

dt_exsitu_period <- lapply(X = datos_time$period_date_strat %>% 
                               unique(),function(x) get_freq_genebank(datos_time, 
                                                                "period_date_strat",
                                                                    x, 
                                                                "cu_variety_name", 
                                                                "data_source",
                                                                c("GRIN database",
                                                                 "CIP's GeneBank"),                                                              "cu_variety_name")
                           ) %>% 
                        rbindlist() %>% as.data.frame()

dt_exsitu_cultivar_period <- lapply(1:nrow(dt_exsitu_period), function(x)  
                                get_exsitu_cultivar(datos_time, dt_exsitu_period[x,], 
                                                    "varie", 
                                                dt_exsitu_period[x,]$period_date_strat)
                                ) %>% rbindlist() %>% as.data.frame()

2.16 Indice de Shanon para analisis de diversidad de los 3 distritos y con las variedades más comunes

library(plyr)
tbl_shanon <-  ddply(tbl_freq_varie_adm3,~ADM3_Name,function(x) {
                        data.frame(SHANNON=diversity(x[-1], index="shannon"))
 })
ggpubr::ggdotchart(tbl_shanon, x = "ADM3_Name", y = "SHANNON", color="ADM3_Name", dot.size = 5, y.text.col=TRUE, rotate=TRUE,  ggtheme = theme_pubr(), palette = c("#00AFBB", "#E7B800", "#FC4E07"), # Custom color palette
                   sorting = "descending" ) + theme_cleveland()  + xlab("Indice de Shannon")

2.17 Indice de Pielou para analisis de equidad de los 3 distritos y con las variedades más comunes

library(plyr)
tbl_pileu_adm3 <-  ddply(tbl_freq_varie_adm3,~ADM3_Name,function(x) {
         data.frame(Pielou=diversity(x[-1], index="simpson")/log(sum(x[-1]>0))) })
flextable(tbl_pileu_adm3)

ADM3_Name

Pielou

Challabamba

0.1727401

Colquepata

0.1739207

Paucartambo

0.1737129

ggpubr::ggdotchart(tbl_pileu_adm3, x = "ADM3_Name", y = "Pielou", color="ADM3_Name", dot.size = 5, y.text.col=TRUE, rotate=TRUE,  ggtheme = theme_pubr(), palette = c("#00AFBB", "#E7B800", "#FC4E07"), # Custom color palette
                   sorting = "descending" ) + theme_cleveland()  + xlab("Indice de Shannon")

2.18 Erosion genética varietal de Colquepata Actulidad vs 2000-2020

out <-  datos_time %>% 
            group_by(ADM3_Name, period_date_strat) %>% 
            dplyr::count(cu_variety_name) %>% 
            relocate("cu_variety_name", .before = "ADM3_Name")

out_colque <- out %>% filter(ADM3_Name=="Colquepata") %>% 
    arrange(desc(n)) %>% 
    pivot_wider(names_from=period_date_strat, values_from= n) %>% 
    relocate(`2021-actualidad`, .before="2001-2020")

out_colque <- out_colque %>% mutate(dif1 = 100- (100*`2021-actualidad`/ `2001-2020`), dif2= 100- (100*`2021-actualidad` / `1981-2000`), dif3= 100- (100*`2021-actualidad` / `1961-1980`), dif4= 100- (100*`2021-actualidad` / `1941-1960`) )

dt_erosion_colque <- out_colque %>% pivot_longer(starts_with("dif"), names_to = "erosion") %>% filter(value>=-100, value<=100) %>%  arrange(desc(value))
#dt_erosion_challabamba[1:10,c(1:4,9)] %>% flextable()

out <- dt_erosion_colque %>% 
            mutate(value=round(value,2)) %>% select(c(1:2,4:5,9))

DT::datatable(out)
#Comparacion actualidad vs 2000-2020
gg_colque <- ggplot(dt_erosion_colque[1:10,] ) + 
    geom_point(aes(x = 1:10, y = `2021-actualidad`),color = "yellow", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = `2001-2020`, ),color = "blue", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = value),color = "grey", size=4,alpha=0.5) +
    #geom_text(aes(label=A), vjust=-1)+
    ylab("Value") +
    xlab("Record Index") +
    theme_minimal()+
    scale_x_continuous(breaks=seq(1,10,1),labels=dt_erosion_colque[1:10,"cu_variety_name"] %>% pull())

gg_colque + ggpubr::rotate_x_text()

2.19 Erosion genética varietal de Colquepata Actulidad vs 1940-1960

 dt_erosion_colque <- out_colque %>% pivot_longer(starts_with("dif"), names_to = "erosion") %>% filter(value>=-100, value<=100) %>% arrange(desc(value))
#dt_erosion_challabamba[1:10,c(1:4,9)] %>% flextable()

out <- dt_erosion_colque %>% 
            mutate(value=round(value,2)) %>% select(c(1:2,4,6,9))

DT::datatable(out)
#Comparacion actualidad vs 2000-2020
gg_colque <- ggplot(dt_erosion_colque[1:10,] ) + 
    geom_point(aes(x = 1:10, y = `2021-actualidad`),color = "yellow", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = `1941-1960`, ),color = "blue", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = value),color = "grey", size=4,alpha=0.5) +
    #geom_text(aes(label=A), vjust=-1)+
    ylab("Value") +
    xlab("Record Index") +
    theme_minimal()+
    scale_x_continuous(breaks=seq(1,10,1),labels=dt_erosion_colque[1:10,"cu_variety_name"] %>% pull())

gg_colque + ggpubr::rotate_x_text()

2.20 Erosion genética varietal de Challabamba Actualidad vs 2000-2020

out <-  datos_time %>% 
            group_by(ADM3_Name, period_date_strat) %>% 
            dplyr::count(cu_variety_name) %>% 
            relocate("cu_variety_name", .before = "ADM3_Name")

out_challabamba <- out %>% filter(ADM3_Name=="Challabamba") %>% 
    arrange(desc(n)) %>% 
    pivot_wider(names_from=period_date_strat, values_from= n) %>% 
    relocate(`2021-actualidad`, .before="2001-2020")

out_challabamba <- out_challabamba %>% mutate(dif1 = 100- (100*`2021-actualidad`/ `2001-2020`), dif2= 100- (100*`2021-actualidad` / `1981-2000`), dif3= 100- (100*`2021-actualidad` / `1961-1980`), dif4= 100- (100*`2021-actualidad` / `1941-1960`) )


dt_erosion_challabamba <- out_challabamba %>% pivot_longer(starts_with("dif"), names_to = "erosion") %>% filter(value>=-100, value<=100) %>% arrange(desc(`2021-actualidad`))
#dt_erosion_challabamba[1:10,c(1:4,9)] %>% flextable()

out <- dt_erosion_challabamba %>% 
            mutate(value=round(value,2)) %>% select(c(1:4,9))

DT::datatable(out)
#Comparacion actualidad vs 2000-2020
gg_challa <- ggplot(dt_erosion_challabamba[1:10,] ) + 
    geom_point(aes(x = 1:10, y = `2021-actualidad`),color = "yellow", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = `2001-2020`, ),color = "blue", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = value),color = "grey", size=4,alpha=0.5) +
    #geom_text(aes(label=A), vjust=-1)+
    ylab("Value") +
    xlab("Record Index") +
    theme_minimal()+
    scale_x_continuous(breaks=seq(1,10,1),labels=dt_erosion_challabamba[1:10,"cu_variety_name"] %>% pull())

gg_challa + ggpubr::rotate_x_text()

2.21 Erosion genética varietal de Paucartambo Actualidad vs 2000-2020

out <-  datos_time %>% 
            group_by(ADM3_Name, period_date_strat) %>% 
            dplyr::count(cu_variety_name) %>% 
            relocate("cu_variety_name", .before = "ADM3_Name")

out_pau <- out %>% filter(ADM3_Name=="Paucartambo") %>% 
    arrange(desc(n)) %>% 
    pivot_wider(names_from=period_date_strat, values_from= n) %>% 
    relocate(`2021-actualidad`, .before="2001-2020")

out_pau <- out_pau %>% mutate(dif1 = 100- (100*`2021-actualidad`/ `2001-2020`), dif2= 100- (100*`2021-actualidad` / `1981-2000`), dif3= 100- (100*`2021-actualidad` / `1961-1980`), dif4= 100- (100*`2021-actualidad` / `1941-1960`) )

dt_erosion_pau <- out_pau %>% pivot_longer(starts_with("dif"), names_to = "erosion") %>% filter(value>=-100, value<=100) %>% arrange(desc(`2021-actualidad`))
#dt_erosion_challabamba[1:10,c(1:4,9)] %>% flextable()

out <- dt_erosion_pau %>% 
            mutate(value=round(value,2)) %>% select(c(1:2,4:5,9))

DT::datatable(out)
#Comparacion actualidad vs 2000-2020
gg_pau <- ggplot(dt_erosion_pau[1:10,] ) + 
    geom_point(aes(x = 1:10, y = `2021-actualidad`),color = "yellow", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = `2001-2020`, ),color = "blue", size=4,alpha=0.5) +
    geom_point(aes(x = 1:10, y = value),color = "grey", size=4,alpha=0.5) +
    #geom_text(aes(label=A), vjust=-1)+
    ylab("Value") +
    xlab("Record Index") +
    theme_minimal()+
    scale_x_continuous(breaks=seq(1,10,1),labels=dt_erosion_pau[1:10,"cu_variety_name"] %>% pull())

gg_pau + ggpubr::rotate_x_text()

2.22 Tendencia temporal por año y categoria varietal

smry_cgroup_year <- janitor::tabyl(datos_time, cu_date_collection, cultivar_group) %>% 
    pivot_longer(2:4, names_to="cultivar_group") %>% 
    group_by(cu_date_collection) %>%  
    mutate(pct=pct_fun(value,2))
#smry_cgroup_year$cu_date_collection <-  lubridate::ym(sprintf("%d-01",smry_cgroup_year$cu_date_collection))

Gráfico de tendencia en el tiempo segun categoria varietal

library(ggthemes)
library(ggHoriPlot)
gho1 <-  smry_cgroup_year %>%  
  ggplot() +
  ggHoriPlot::geom_horizon(aes(cu_date_collection, 
                   pct,fill = ..Cutpoints..), 
                   origin = 0, horizonscale = c(0,20,40,60,80,100)) +
          scale_fill_hcl(palette = 'RdBu', reverse = T) +
        facet_grid(cultivar_group~.) +
  theme_few() +
  theme(
    panel.spacing.y=unit(0, "lines"),
    strip.text.y = element_text(angle = 0, hjust = 0),
    legend.position = 'top',
    #axis.text.y = element_blank(),
    #axis.title.y = element_blank(),
    #axis.ticks.y = element_blank(),
    panel.border = element_blank()
    )  +
    scale_x_continuous(breaks = smry_cgroup_year$cu_date_collection, labels = smry_cgroup_year$cu_date_collection) + ggpubr::rotate_x_text()
gho1    

2.22.1 Tendencia temporal por estratificacion #2 y categoria varietal

smry_cgroup_period <- janitor::tabyl(datos_time, period_date_strat, cultivar_group) %>%
        pivot_longer(2:4, names_to="cultivar_group") %>% 
        group_by(period_date_strat) %>%  
        mutate(pct=pct_fun(value,2))
library(ggthemes)
library(ggHoriPlot)
gho2 <- smry_cgroup_period %>%  
  ggplot() +
  ggHoriPlot::geom_horizon(aes(period_date_strat, 
                   pct,fill = ..Cutpoints..), 
                   origin = 0, horizonscale = c(0,20,40,60,80,100)) +
  scale_fill_hcl(palette = 'RdBu', reverse = T) +
  facet_grid(cultivar_group~.) +
  theme_few() +
  theme(
    panel.spacing.y=unit(0, "lines"),
    strip.text.y = element_text(angle = 0, hjust = 0),
    legend.position = 'top',
    #axis.text.y = element_blank(),
    #axis.title.y = element_blank(),
    #axis.ticks.y = element_blank(),
    panel.border = element_blank()
    )  +
    ggpubr::rotate_x_text()
gho2

2.23 Relacion entre zonas de altitud vs periodo de años y cultivar_group

a2 <- ggplot(datos_time, aes(x = period_date_strat, y = elevation, color=cultivar_group)) + geom_boxplot() + theme_bw() + stat_summary(fun.y=mean, geom="point", shape=20, size=1, color="yellow", fill="yellow") + facet_wrap(~ADM3_Name)
a2 + ggpubr::rotate_x_text()

lm1 <- lm(elevation~ period_date_strat + cultivar_group, data= datos_time)
aov_lm1 <- aov(lm1)
summary(aov_lm1)
                     Df    Sum Sq  Mean Sq F value   Pr(>F)    
period_date_strat     4 8.028e+07 20068964  231.79  < 2e-16 ***
cultivar_group        2 2.224e+06  1111773   12.84 2.68e-06 ***
Residuals         15250 1.320e+09    86581                     
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
out <- HSD.test(aov_lm1,"period_date_strat", group=TRUE,console=TRUE)

Study: aov_lm1 ~ "period_date_strat"

HSD Test for elevation 

Mean Square Error:  86580.61 

period_date_strat,  means

                elevation      std    r     Min     Max
1941-1960        3598.765 288.1917  434 2721.68 3917.47
1961-1980        3312.494 671.9209  308 2233.06 4254.93
1981-2000        3649.146 405.7562 1243 2721.68 4425.53
2001-2020        3715.493 201.7777 5356 2721.68 4138.10
2021-actualidad  3766.012 302.7576 7916 2782.43 4518.66

Alpha: 0.05 ; DF Error: 15250 
Critical Value of Studentized Range: 3.858118 

Groups according to probability of means differences and alpha level( 0.05 )

Treatments with the same letter are not significantly different.

                elevation groups
2021-actualidad  3766.012      a
2001-2020        3715.493      b
1981-2000        3649.146      c
1941-1960        3598.765      d
1961-1980        3312.494      e
DT::datatable(out$means)

Study: aov_lm1 ~ "cultivar_group"

HSD Test for elevation 

Mean Square Error:  86580.61 

cultivar_group,  means

                elevation      std     r     Min     Max
Bitter landrace  3752.832 294.4079  1994 2233.06 4518.66
Floury landrace  3721.900 309.0494 12155 2233.06 4518.66
Modern variety   3706.769 246.4412  1108 2782.43 4473.75

Alpha: 0.05 ; DF Error: 15250 
Critical Value of Studentized Range: 3.314818 

Groups according to probability of means differences and alpha level( 0.05 )

Treatments with the same letter are not significantly different.

                elevation groups
Bitter landrace  3752.832      a
Floury landrace  3721.900      b
Modern variety   3706.769      b