#1.LIBRARIES (PACKAGES "R") ===========================================================================================================================
  library(openxlsx)
  library(reshape2)
  library(stringr)
  library(plyr)

#2.PROCESSING DATABASE FILE "CROPS" - NO_HUMIDITY - EUROSTAT ========================================================================================== 
  #Data import and transpose data by year - Variable AREA
  ELENCO_sheets <- c("Sheet 3", "Sheet 4", "Sheet 5", "Sheet 6", "Sheet 7", "Sheet 8", "Sheet 9",
                     "Sheet 10", "Sheet 11", "Sheet 12", "Sheet 13", "Sheet 14", "Sheet 15", "Sheet 16", "Sheet 17", 
                     "Sheet 18", "Sheet 19", "Sheet 20", "Sheet 21", "Sheet 22", "Sheet 23", "Sheet 24", "Sheet 25", 
                     "Sheet 26", "Sheet 27", "Sheet 28", "Sheet 29", "Sheet 30", "Sheet 31", "Sheet 32", "Sheet 33", 
                     "Sheet 34", "Sheet 35", "Sheet 36", "Sheet 37", "Sheet 38", "Sheet 39", "Sheet 40", "Sheet 41", 
                     "Sheet 42", "Sheet 43", "Sheet 44", "Sheet 45", "Sheet 46", "Sheet 47", "Sheet 48", "Sheet 49",
                     "Sheet 50", "Sheet 51", "Sheet 52", "Sheet 53", "Sheet 54", "Sheet 55", "Sheet 56", "Sheet 57",
                     "Sheet 58", "Sheet 59", "Sheet 60", "Sheet 61", "Sheet 62", "Sheet 63", "Sheet 64", "Sheet 65",
                     "Sheet 66", "Sheet 67", "Sheet 68", "Sheet 69", "Sheet 70", "Sheet 71", "Sheet 72", "Sheet 73",
                     "Sheet 74", "Sheet 75", "Sheet 76", "Sheet 77", "Sheet 78", "Sheet 79")
  length(ELENCO_sheets)
  ELENCO_colonne <- seq(3,19,2)
  
  setwd("C:\\Users\\................\\data folder")
  DB_SUPERFICI_final <- data.frame()
  
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 77 fogli", sep = ""))
       n_sheet <- ELENCO_sheets[i]
       coltura <- colnames(read.xlsx(xlsxFile = "SUP_apro_cpshr__custom_4606765_spreadsheet.xlsx", sheet = n_sheet, rows = 6, cols = 3))  
       coltura <- gsub(".", " ", coltura, fixed = TRUE)
       table <- read.xlsx(xlsxFile = "SUP_apro_cpshr__custom_4606765_spreadsheet.xlsx", sheet = n_sheet, 
                          rows = c(10:647), cols = c(1:20))
       colnames(table)[1:20] <- c("GEO_CODE", "GEO_TABLES", "a_2013", "flag_13", "a_2014", "flag_14", "a_2015", "flag_15",
                                  "a_2016", "flag_16", "a_2017", "flag_17", "a_2018", "flag_18", "a_2019", "flag_19",
                                  "a_2020", "flag_20", "a_2021", "flag_21")
       
       #2013
       temp_1a <- table[,c(1:3)]
       temp_1a <- melt(temp_1a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1b <- table[,c(1,2,4)]
       temp_1b <- melt(temp_1b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1 <- data.frame(temp_1a, flag = temp_1b[,4])
       remove(temp_1a, temp_1b)
       #2014
       temp_2a <- table[,c(1,2,5)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,6)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,7)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,8)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,9)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,10)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,11)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,12)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
       remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,13)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,14)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,15)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,16)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,17)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,18)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
       remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,19)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,20)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
       
       TEMP <- rbind(temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       TEMP$Coltura <- coltura
       DB_SUPERFICI_final <- rbind(DB_SUPERFICI_final, TEMP)}

       remove(i, table, coltura, n_sheet, TEMP, temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       remove(ELENCO_sheets, ELENCO_colonne)    
 
  colnames(DB_SUPERFICI_final)[3:6] <- c("ANNO", "SUPERFICIE", "FLAG_SUP", "COLTURA")
  DB_SUPERFICI_final$ANNO <- gsub("a_", "", DB_SUPERFICI_final$ANNO)
  DB_SUPERFICI_final <- DB_SUPERFICI_final[,c(3,1,2,6,4,5)]
  DB_SUPERFICI_final$SUPERFICIE_NUM <- as.numeric(DB_SUPERFICI_final$SUPERFICIE)
  
  DB_SUPERFICI_final$CODICE_COLTURA <- ""
  for (i in 1:nrow(DB_SUPERFICI_final)){
       word.list <- str_split(DB_SUPERFICI_final$COLTURA[i], "\\s+")
       words     <- unlist(word.list)
       cod_colt  <- words[length(words)] 
       DB_SUPERFICI_final$CODICE_COLTURA[i] <- cod_colt
       print(i)}
       remove(i, word.list, words, cod_colt)
  
  #Data import and transpose data by year - Variable PRODUCTION      
  ELENCO_sheets <- c("Sheet 3", "Sheet 4", "Sheet 5", "Sheet 6", "Sheet 7", "Sheet 8", "Sheet 9",
                     "Sheet 10", "Sheet 11", "Sheet 12", "Sheet 13", "Sheet 14", "Sheet 15", "Sheet 16", "Sheet 17", 
                     "Sheet 18", "Sheet 19", "Sheet 20", "Sheet 21", "Sheet 22", "Sheet 23", "Sheet 24", "Sheet 25", 
                     "Sheet 26", "Sheet 27", "Sheet 28", "Sheet 29", "Sheet 30", "Sheet 31", "Sheet 32", "Sheet 33", 
                     "Sheet 34", "Sheet 35", "Sheet 36", "Sheet 37", "Sheet 38", "Sheet 39", "Sheet 40", "Sheet 41", 
                     "Sheet 42", "Sheet 43", "Sheet 44", "Sheet 45", "Sheet 46", "Sheet 47", "Sheet 48", "Sheet 49",
                     "Sheet 50", "Sheet 51", "Sheet 52", "Sheet 53", "Sheet 54", "Sheet 55", "Sheet 56", "Sheet 57",
                     "Sheet 58", "Sheet 59", "Sheet 60", "Sheet 61", "Sheet 62", "Sheet 63", "Sheet 64", "Sheet 65",
                     "Sheet 66", "Sheet 67", "Sheet 68", "Sheet 69", "Sheet 70", "Sheet 71", "Sheet 72", "Sheet 73",
                     "Sheet 74", "Sheet 75", "Sheet 76", "Sheet 77", "Sheet 78", "Sheet 79")
  length(ELENCO_sheets)
  ELENCO_colonne <- seq(3,19,2)
   
  setwd("C:\\Users\\\\................\\data folder")
  DB_PRODUZIONI_final <- data.frame()
       
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 77 fogli", sep = ""))
       n_sheet <- ELENCO_sheets[i]
       coltura <- colnames(read.xlsx(xlsxFile = "PROD_apro_cpshr__custom_4606814_spreadsheet.xlsx", sheet = n_sheet, rows = 6, cols = 3))  
       coltura <- gsub(".", " ", coltura, fixed = TRUE)
       table <- read.xlsx(xlsxFile = "PROD_apro_cpshr__custom_4606814_spreadsheet.xlsx", sheet = n_sheet, 
                          rows = c(10:647), cols = c(1:20))
       colnames(table)[1:20] <- c("GEO_CODE", "GEO_TABLES", "a_2013", "flag_13", "a_2014", "flag_14", "a_2015", "flag_15",
                                  "a_2016", "flag_16", "a_2017", "flag_17", "a_2018", "flag_18", "a_2019", "flag_19",
                                  "a_2020", "flag_20", "a_2021", "flag_21")
       
       #2013
       temp_1a <- table[,c(1:3)]
       temp_1a <- melt(temp_1a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1b <- table[,c(1,2,4)]
       temp_1b <- melt(temp_1b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1 <- data.frame(temp_1a, flag = temp_1b[,4])
       remove(temp_1a, temp_1b)
       #2014
       temp_2a <- table[,c(1,2,5)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,6)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,7)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,8)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,9)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,10)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,11)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,12)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
       remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,13)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,14)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,15)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,16)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,17)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,18)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
       remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,19)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,20)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
       
       TEMP <- rbind(temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       TEMP$Coltura <- coltura
       DB_PRODUZIONI_final <- rbind(DB_PRODUZIONI_final, TEMP)}
  
  remove(i, table, coltura, n_sheet, TEMP, temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  remove(ELENCO_sheets, ELENCO_colonne)
   
  colnames(DB_PRODUZIONI_final)[3:6] <- c("ANNO", "PRODUZIONE", "FLAG_PROD", "COLTURA")
  DB_PRODUZIONI_final$ANNO <- gsub("a_", "", DB_PRODUZIONI_final$ANNO)
  DB_PRODUZIONI_final <- DB_PRODUZIONI_final[,c(3,1,2,6,4,5)]
  DB_PRODUZIONI_final$PRODUZIONE_NUM <- as.numeric(DB_PRODUZIONI_final$PRODUZIONE)
  
  DB_PRODUZIONI_final$CODICE_COLTURA <- ""
  for (i in 1:nrow(DB_PRODUZIONI_final)){
       word.list <- str_split(DB_PRODUZIONI_final$COLTURA[i], "\\s+")
       words     <- unlist(word.list)
       cod_colt  <- words[length(words)] 
       DB_PRODUZIONI_final$CODICE_COLTURA[i] <- cod_colt
       print(i)}
       remove(i, word.list, words, cod_colt)

  #Data import and transpose data by year - Variable UTILISED AGRICULTURAL AREA (UAA)     
  ELENCO_sheets <- c("Sheet 1")
      
  setwd("C:\\Users\\................\\data folder")
  
  DB_SAU_final <- data.frame()
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 77 fogli", sep = ""))
       n_sheet <- ELENCO_sheets[i]
       coltura <- colnames(read.xlsx(xlsxFile = "SAU_apro_cpshr__custom_4606845_spreadsheet.xlsx", sheet = n_sheet, rows = 6, cols = 3))  
       coltura <- gsub(".", " ", coltura, fixed = TRUE)
       table <- read.xlsx(xlsxFile = "SAU_apro_cpshr__custom_4606845_spreadsheet.xlsx", sheet = n_sheet, 
                          rows = c(10:647), cols = c(1:20))
       colnames(table)[1:20] <- c("GEO_CODE", "GEO_TABLES", "a_2013", "flag_13", "a_2014", "flag_14", "a_2015", "flag_15",
                                  "a_2016", "flag_16", "a_2017", "flag_17", "a_2018", "flag_18", "a_2019", "flag_19",
                                  "a_2020", "flag_20", "a_2021", "flag_21")
     
       #2013
       temp_1a <- table[,c(1:3)]
       temp_1a <- melt(temp_1a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1b <- table[,c(1,2,4)]
       temp_1b <- melt(temp_1b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1 <- data.frame(temp_1a, flag = temp_1b[,4])
       remove(temp_1a, temp_1b)
       #2014
       temp_2a <- table[,c(1,2,5)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,6)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,7)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,8)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,9)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,10)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,11)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,12)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
      remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,13)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,14)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,15)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,16)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,17)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,18)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
     remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,19)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,20)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
     
  TEMP <- rbind(temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  TEMP$Coltura <- coltura
  DB_SAU_final <- rbind(DB_SAU_final, TEMP)}
  remove(i, table, coltura, n_sheet, TEMP, temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  remove(ELENCO_sheets)
     
  colnames(DB_SAU_final)[3:6] <- c("ANNO", "SAU", "FLAG_SAU", "COLTURA")
  DB_SAU_final$ANNO <- gsub("a_", "", DB_SAU_final$ANNO)
  DB_SAU_final <- DB_SAU_final[,c(3,1,2,6,4,5)]
  DB_SAU_final$SAU_NUM <- as.numeric(DB_SAU_final$SAU)
  
  DB_SAU_final$CODICE_COLTURA <- ""
  for (i in 1:nrow(DB_SAU_final)){
       word.list <- str_split(DB_SAU_final$COLTURA[i], "\\s+")
       words     <- unlist(word.list)
       cod_colt  <- words[length(words)] 
       DB_SAU_final$CODICE_COLTURA[i] <- cod_colt
       print(i)}
  remove(i, word.list, words, cod_colt)
  
  #Aggregation Database and information alignment
  DATABASE_1 <- merge(DB_SUPERFICI_final, DB_PRODUZIONI_final[,c(1,2,8,5:7)], by = c("ANNO","GEO_CODE","CODICE_COLTURA"), all.x = TRUE, all.y = TRUE)
  DATABASE_1 <- merge(DATABASE_1, DB_SAU_final[,c(1,2,5:7)], by = c("ANNO","GEO_CODE"), all.x = TRUE, all.y = FALSE)
  DATABASE_1 <- subset(DATABASE_1, select = c("ANNO","GEO_CODE","GEO_TABLES",
                                              "CODICE_COLTURA","COLTURA","SUPERFICIE_NUM","FLAG_SUP",
                                              "PRODUZIONE_NUM","FLAG_PROD","SAU_NUM","FLAG_SAU"))
  colnames(DATABASE_1)[1:11] <- c("YEAR","NUTS2_CODE","NUTS2_NAME",
                                  "CROP_CODE","CROP_NAME","AREA","FLAG_AREA",
                                  "PROD","FLAG_PROD","SAU","FLAG_SAU")

  DATABASE_1$SOURCE_CODE <- "[APRO_CPSHR]"
  DATABASE_1$SOURCE_LINK <- "https://ec.europa.eu/eurostat/databrowser/view/apro_cpshr/default/table?lang=en"
  DATABASE_1$SOURCE_NAME <- "Crop production in EU standard humidity by NUTS 2 regions"
  
  #Selection of EU countries and reference period (year 2017-2021)
  DATABASE_1$NUTS0_CODE <- substr(DATABASE_1$NUTS2_CODE,1,2)
  CODE_EU_27 <- c("BE",	"BG",	"CZ",	"DK",	"DE",	"EE",	"IE",	"EL",	"ES",	"FR",	"HR",	"IT",	
                  "CY",	"LV",	"LT",	"LU",	"HU",	"MT",	"NL",	"AT",	"PL",	"PT",	"RO",	"SI",	
                  "SK",	"FI",	"SE")
  DATABASE_1 <- DATABASE_1[DATABASE_1$NUTS0_CODE %in% CODE_EU_27,]
  remove(CODE_EU_27) 
  DATABASE_1 <- subset(DATABASE_1, DATABASE_1$YEAR>=2017)
  
  temp_TAB <- unique(DATABASE_1[,c(2,3)])
  temp_TAB <- subset(temp_TAB, nchar(temp_TAB$NUTS2_CODE)==2)
  colnames(temp_TAB)[1:2] <- c("NUTS0_CODE","NUTS0_NAME")
  DATABASE_1 <- merge(DATABASE_1, temp_TAB, by = "NUTS0_CODE", all.x = TRUE, all.y = FALSE)
  remove(temp_TAB)
  
  DATABASE_1$CROP_CODE <- gsub("\\[|\\]","", DATABASE_1$CROP_CODE)
  
  remove(DB_SUPERFICI_final, DB_PRODUZIONI_final, DB_SAU_final)
  
#3.PROCESSING DATABASE FILE "CROPS" - WITH_HUMIDITY - EUROSTAT =======================================================================================
  #Data import and transpose data by year - Variable AREA
  ELENCO_sheets <- getSheetNames("C:\\Users\\...............\\AREA_apro_cpnh1__custom_5395178_spreadsheet.xlsx")
  ELENCO_sheets <- ELENCO_sheets[5:199]
  ELENCO_sheets
  
  ELENCO_colonne <- seq(3,19,2)
  setwd("C:\\Users\\.................\\data folder")
  
  DB_SUPERFICI_final <- data.frame()
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 197 fogli", sep = ""))
       n_sheet <- ELENCO_sheets[i]
       coltura <- colnames(read.xlsx(xlsxFile = "AREA_apro_cpnh1__custom_5395178_spreadsheet.xlsx", sheet = n_sheet, rows = 6, cols = 3))  
       coltura <- gsub(".", " ", coltura, fixed = TRUE)
       table <- read.xlsx(xlsxFile = "AREA_apro_cpnh1__custom_5395178_spreadsheet.xlsx", sheet = n_sheet, 
                          rows = c(10:51), cols = c(1:20))
       colnames(table)[1:20] <- c("GEO_CODE", "GEO_TABLES", "a_2013", "flag_13", "a_2014", "flag_14", "a_2015", 
                                  "flag_15","a_2016", "flag_16", "a_2017", "flag_17", "a_2018", "flag_18", 
                                  "a_2019", "flag_19","a_2020", "flag_20", "a_2021", "flag_21")
       
       #2013
       temp_1a <- table[,c(1,2,3)]
       temp_1a <- melt(temp_1a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1b <- table[,c(1,2,4)]
       temp_1b <- melt(temp_1b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1 <- data.frame(temp_1a, flag = temp_1b[,4])
       remove(temp_1a, temp_1b)
       #2014
       temp_2a <- table[,c(1,2,5)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,6)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,7)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,8)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,9)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,10)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,11)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,12)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
       remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,13)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,14)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,15)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,16)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,17)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,18)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
       remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,19)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,20)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
       
       TEMP <- rbind(temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       TEMP$Coltura <- coltura
       DB_SUPERFICI_final <- rbind(DB_SUPERFICI_final, TEMP)}
  
  remove(i, table, coltura, n_sheet, TEMP, temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  remove(ELENCO_sheets, ELENCO_colonne)
  
  colnames(DB_SUPERFICI_final)[3:6] <- c("ANNO", "AREA", "FLAG_AREA", "CROP_NAME")
  DB_SUPERFICI_final$ANNO <- gsub("a_", "", DB_SUPERFICI_final$ANNO)
  DB_SUPERFICI_final <- DB_SUPERFICI_final[,c(3,1,2,6,4,5)]
  DB_SUPERFICI_final$AREA_NUM <- as.numeric(DB_SUPERFICI_final$AREA)
  
  DB_SUPERFICI_final$CROP_CODE <- ""
  for (i in 1:nrow(DB_SUPERFICI_final)){
       word.list <- str_split(DB_SUPERFICI_final$CROP_NAME[i], "\\s+")
       words     <- unlist(word.list)
       cod_colt  <- words[length(words)] 
       DB_SUPERFICI_final$CROP_CODE[i] <- cod_colt
       print(i)}
       remove(i, word.list, words, cod_colt)
  
  DB_SUPERFICI_final <- DB_SUPERFICI_final[,c(1:3,8,4,7,6)]
  colnames(DB_SUPERFICI_final)[c(2,3,6,7)] <- c("NUTS0_CODE","NUTS0_NAME","AREA","FLAG_AREA")
  DB_SUPERFICI_final <- DB_SUPERFICI_final[order(DB_SUPERFICI_final$ANNO, DB_SUPERFICI_final$NUTS0_CODE, DB_SUPERFICI_final$CROP_CODE),]
  
  DB_SUPERFICI_final <- subset(DB_SUPERFICI_final, DB_SUPERFICI_final$ANNO>=2017)
  
  #Data import and transpose data by year - Variable PRODUCTION      
  ELENCO_sheets <- getSheetNames("C:\\Users\\........\\PROD_apro_cpnh1__custom_5396069_spreadsheet.xlsx")
  ELENCO_sheets <- ELENCO_sheets[5:199]
  ELENCO_sheets
  
  length(ELENCO_sheets)
  ELENCO_colonne <- seq(3,19,2)
  
  setwd("C:\\Users\\................\\data folder")
  
  DB_PRODUZIONI_final <- data.frame()
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 197 fogli", sep = ""))
       n_sheet <- ELENCO_sheets[i]
       coltura <- colnames(read.xlsx(xlsxFile = "PROD_apro_cpnh1__custom_5396069_spreadsheet.xlsx", sheet = n_sheet, rows = 6, cols = 3))  
       coltura <- gsub(".", " ", coltura, fixed = TRUE)
       table <- read.xlsx(xlsxFile = "PROD_apro_cpnh1__custom_5396069_spreadsheet.xlsx", sheet = n_sheet, 
                          rows = c(10:51), cols = c(1:20))
       colnames(table)[1:20] <- c("GEO_CODE", "GEO_TABLES", "a_2013", "flag_13", "a_2014", "flag_14", "a_2015", "flag_15",
                                  "a_2016", "flag_16", "a_2017", "flag_17", "a_2018", "flag_18", "a_2019", "flag_19",
                                  "a_2020", "flag_20", "a_2021", "flag_21")
       
       #2013
       temp_1a <- table[,c(1:3)]
       temp_1a <- melt(temp_1a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1b <- table[,c(1,2,4)]
       temp_1b <- melt(temp_1b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_1 <- data.frame(temp_1a, flag = temp_1b[,4])
       remove(temp_1a, temp_1b)
       #2014
       temp_2a <- table[,c(1,2,5)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,6)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,7)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,8)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,9)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,10)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,11)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,12)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
       remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,13)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,14)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,15)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,16)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,17)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,18)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
       remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,19)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,20)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
       
       TEMP <- rbind(temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       TEMP$Coltura <- coltura
       DB_PRODUZIONI_final <- rbind(DB_PRODUZIONI_final, TEMP)}
  
  remove(i, table, coltura, n_sheet, TEMP, temp_1, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  remove(ELENCO_sheets, ELENCO_colonne)
  
  colnames(DB_PRODUZIONI_final)[3:6] <- c("ANNO", "PRODUCTION", "FLAG_PROD", "CROP_NAME")
  DB_PRODUZIONI_final$ANNO <- gsub("a_", "", DB_PRODUZIONI_final$ANNO)
  DB_PRODUZIONI_final <- DB_PRODUZIONI_final[,c(3,1,2,6,4,5)]
  DB_PRODUZIONI_final$PRODUCTION_NUM <- as.numeric(DB_PRODUZIONI_final$PRODUCTION)
  
  DB_PRODUZIONI_final$CROP_CODE <- ""
  for (i in 1:nrow(DB_PRODUZIONI_final)){
       word.list <- str_split(DB_PRODUZIONI_final$CROP_NAME[i], "\\s+")
       words     <- unlist(word.list)
       cod_colt  <- words[length(words)] 
       DB_PRODUZIONI_final$CROP_CODE[i] <- cod_colt
       print(i)}
  remove(i, word.list, words, cod_colt)
  
  DB_PRODUZIONI_final <- DB_PRODUZIONI_final[,c(1:3,8,4,7,6)]
  colnames(DB_PRODUZIONI_final)[c(2,3,6,7)] <- c("NUTS0_CODE","NUTS0_NAME","PRODUCTION","FLAG_PROD")
  DB_PRODUZIONI_final <- DB_PRODUZIONI_final[order(DB_PRODUZIONI_final$ANNO, DB_PRODUZIONI_final$NUTS0_CODE, DB_PRODUZIONI_final$CROP_CODE),]
  
  DB_PRODUZIONI_final <- subset(DB_PRODUZIONI_final, DB_PRODUZIONI_final$ANNO>=2017)
  
  #Data import and transpose data by year - Variable UTILISED AGRICULTURAL AREA (UAA)     
  ELENCO_sheets <- c("Sheet 1")
  
  setwd("C:\\Users\\................\\data folder")
  
  DB_SAU_final <- data.frame()
  for (i in 1:length(ELENCO_sheets)){
       print(paste(ELENCO_sheets[i], " - Caricamento del foglio ", i, " di 1 foglio", sep = ""))
       table <- read.xlsx(xlsxFile = "SAU_apro_cpnh1__custom_5390220_spreadsheet.xlsx", sheet = ELENCO_sheets, 
                          rows = c(10:51), cols = c(1:18))
       colnames(table)[1:18] <- c("GEO_CODE", "GEO_TABLES", "a_2014", "flag_14", "a_2015", "flag_15","a_2016", "flag_16",
                                  "a_2017", "flag_17", "a_2018", "flag_18", "a_2019", "flag_19","a_2020", "flag_20",
                                  "a_2021", "flag_21")
       
       #2014
       temp_2a <- table[,c(1,2,3)]
       temp_2a <- melt(temp_2a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2b <- table[,c(1,2,4)]
       temp_2b <- melt(temp_2b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_2 <- data.frame(temp_2a, flag = temp_2b[,4])
       remove(temp_2a, temp_2b)
       #2015
       temp_3a <- table[,c(1,2,5)]
       temp_3a <- melt(temp_3a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3b <- table[,c(1,2,6)]
       temp_3b <- melt(temp_3b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_3 <- data.frame(temp_3a, flag = temp_3b[,4])
       remove(temp_3a, temp_3b)
       #2016
       temp_4a <- table[,c(1,2,7)]
       temp_4a <- melt(temp_4a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4b <- table[,c(1,2,8)]
       temp_4b <- melt(temp_4b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_4 <- data.frame(temp_4a, flag = temp_4b[,4])
       remove(temp_4a, temp_4b)
       #2017
       temp_5a <- table[,c(1,2,9)]
       temp_5a <- melt(temp_5a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5b <- table[,c(1,2,10)]
       temp_5b <- melt(temp_5b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_5 <- data.frame(temp_5a, flag = temp_5b[,4])
       remove(temp_5a, temp_5b)
       #2018
       temp_6a <- table[,c(1,2,11)]
       temp_6a <- melt(temp_6a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6b <- table[,c(1,2,12)]
       temp_6b <- melt(temp_6b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_6 <- data.frame(temp_6a, flag = temp_6b[,4])
       remove(temp_6a, temp_6b)
       #2019
       temp_7a <- table[,c(1,2,13)]
       temp_7a <- melt(temp_7a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7b <- table[,c(1,2,14)]
       temp_7b <- melt(temp_7b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_7 <- data.frame(temp_7a, flag = temp_7b[,4])
       remove(temp_7a, temp_7b)
       #2020
       temp_8a <- table[,c(1,2,15)]
       temp_8a <- melt(temp_8a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8b <- table[,c(1,2,16)]
       temp_8b <- melt(temp_8b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_8 <- data.frame(temp_8a, flag = temp_8b[,4])
       remove(temp_8a, temp_8b)
       #2021
       temp_9a <- table[,c(1,2,17)]
       temp_9a <- melt(temp_9a, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9b <- table[,c(1,2,18)]
       temp_9b <- melt(temp_9b, id = c("GEO_CODE", "GEO_TABLES"))
       temp_9 <- data.frame(temp_9a, flag = temp_9b[,4])
       remove(temp_9a, temp_9b)
       
       TEMP <- rbind(temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
       DB_SAU_final <- rbind(DB_SAU_final, TEMP)}
     
  remove(i, table, TEMP, temp_2, temp_3, temp_4, temp_5, temp_6, temp_7, temp_8, temp_9)
  remove(ELENCO_sheets)
  
  colnames(DB_SAU_final)[3:5] <- c("ANNO", "UAA", "FLAG_UAA")
  DB_SAU_final$ANNO <- gsub("a_", "", DB_SAU_final$ANNO)
  DB_SAU_final <- DB_SAU_final[,c(3,1,2,4,5)]
  DB_SAU_final$UAA_NUM <- as.numeric(DB_SAU_final$UAA)
  
  DB_SAU_final <- DB_SAU_final[,c(1:3,6,5)]
  colnames(DB_SAU_final)[c(2,3,4,5)] <- c("NUTS0_CODE","NUTS0_NAME","UAA","FLAG_UAA")
  DB_SAU_final <- DB_SAU_final[order(DB_SAU_final$ANNO, DB_SAU_final$NUTS0_CODE),]
  
  DB_SAU_final <- subset(DB_SAU_final, DB_SAU_final$ANNO>=2017)    

  #Aggregation Database and information alignments
  DB_SUPERFICI_final$AREA        <- as.numeric(DB_SUPERFICI_final$AREA)
  DB_PRODUZIONI_final$PRODUCTION <- as.numeric(DB_PRODUZIONI_final$PRODUCTION)
  DB_SAU_final$UAA               <- as.numeric(DB_SAU_final$UAA)
  
  DATABASE_2 <- merge(DB_SUPERFICI_final, DB_PRODUZIONI_final[,c(1,2,4,6,7)], by = c("ANNO","NUTS0_CODE","CROP_CODE"), all.x = TRUE, all.y = TRUE)
  DATABASE_2 <- merge(DATABASE_2, DB_SAU_final[,c(1,2,4,5)], by = c("ANNO","NUTS0_CODE"), all.x = TRUE, all.y = FALSE)
  colnames(DATABASE_2)
  DATABASE_2 <- subset(DATABASE_2, select = c("ANNO","NUTS0_CODE","NUTS0_NAME",
                                              "CROP_CODE","CROP_NAME","AREA","FLAG_AREA",
                                              "PRODUCTION","FLAG_PROD","UAA","FLAG_UAA"))
  
  colnames(DATABASE_2)
  colnames(DATABASE_2)[1:11] <- c("YEAR","NUTS0_CODE","NUTS0_NAME",
                                  "CROP_CODE","CROP_NAME","AREA","FLAG_AREA",
                                  "PROD","FLAG_PROD","SAU","FLAG_SAU")

  DATABASE_2$SOURCE_CODE <- "[APRO_CPNH1]"
  DATABASE_2$SOURCE_LINK <- "https://ec.europa.eu/eurostat/databrowser/view/apro_cpnh1/default/table?lang=en"
  DATABASE_2$SOURCE_NAME <- "Crop production in national humidity"
  
  DATABASE_2$NUTS2_CODE <- DATABASE_2$NUTS0_CODE
  DATABASE_2$NUTS2_NAME <- DATABASE_2$NUTS0_NAME
  
  DATABASE_2 <- subset(DATABASE_2, select = c("YEAR","NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME",
                                              "CROP_CODE","CROP_NAME","AREA","FLAG_AREA", "PROD","FLAG_PROD",
                                              "SAU","FLAG_SAU","SOURCE_CODE","SOURCE_LINK","SOURCE_NAME"))
  
  DATABASE_2$CROP_CODE <- gsub("\\[|\\]","", DATABASE_2$CROP_CODE)
  
  CODE_EU_27 <- c("BE",	"BG",	"CZ",	"DK",	"DE",	"EE",	"IE",	"EL",	"ES",	"FR",	"HR",	"IT",	
                  "CY",	"LV",	"LT",	"LU",	"HU",	"MT",	"NL",	"AT",	"PL",	"PT",	"RO",	"SI",	
                  "SK",	"FI",	"SE")
  DATABASE_2 <- DATABASE_2[DATABASE_2$NUTS0_CODE %in% CODE_EU_27,]
  remove(CODE_EU_27) 
  
  remove(DB_PRODUZIONI_final, DB_SUPERFICI_final, DB_SAU_final)
  
  #Export temporary Database_1 and Datanase_2 created
  setwd("C:\\Users\\................\\data folder")
  wb <- createWorkbook()
  addWorksheet(wb, sheetName = "DATABASE_1") 
  writeData(wb, sheet = 1, x = DATABASE_1, na.string = "")        #keep.na = TRUE
  addWorksheet(wb, sheetName = "DATABASE_2") 
  writeData(wb, sheet = 2, x = DATABASE_2, na.string = "")
  saveWorkbook(wb, "temp_DATABASE_1_2.xlsx", overwrite = TRUE)
  remove(wb)
  
#4.APPENDING DATABASE_1 (NO_HUMIDITY) & DATABASE_2 (WITH_HUMIDITY) ==================================================================================
  rm(list=ls())
  
  library(openxlsx)
  library(reshape2)
  library(plyr)
  
  DATABASE_1 <- read.xlsx("C:\\Users\\......\\temp_DATABASE_1_2.xlsx", sheet = "DATABASE_1")
  DATABASE_2 <- read.xlsx("C:\\Users\\......\\temp_DATABASE_1_2.xlsx", sheet = "DATABASE_2")
  
  DATABASE_1 <- subset(DATABASE_1, select = c("YEAR","NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME",
                                              "CROP_CODE","CROP_NAME","AREA","FLAG_AREA", "PROD","FLAG_PROD",
                                              "SAU","FLAG_SAU","SOURCE_CODE","SOURCE_LINK","SOURCE_NAME"))
  DATABASE_2 <- subset(DATABASE_2, select = c("YEAR","NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME",
                                              "CROP_CODE","CROP_NAME","AREA","FLAG_AREA", "PROD","FLAG_PROD",
                                              "SAU","FLAG_SAU","SOURCE_CODE","SOURCE_LINK","SOURCE_NAME"))
  
  TAB_1 <- data.frame(YEAR = DATABASE_1$YEAR, NUTS2_CODE = DATABASE_1$NUTS2_CODE, CROP_CODE = DATABASE_1$CROP_CODE, FLAG_noHUM = 1)
  TAB_2 <- data.frame(YEAR = DATABASE_2$YEAR, NUTS2_CODE = DATABASE_2$NUTS2_CODE, CROP_CODE = DATABASE_2$CROP_CODE, FLAG_HUM = 1)
  
  temp_TAB <- merge(TAB_1, TAB_2, by = c("YEAR","NUTS2_CODE","CROP_CODE"), all.x = TRUE, all.y = TRUE)
  remove(TAB_1, TAB_2)
  temp_TAB$KEEP_ROW <- 0
  temp_TAB$KEEP_ROW <- ifelse(temp_TAB$FLAG_noHUM==1, 1, temp_TAB$KEEP_ROW) 
  temp_TAB$KEEP_ROW <- ifelse(is.na(temp_TAB$FLAG_noHUM)==TRUE & temp_TAB$FLAG_HUM==1, 2, temp_TAB$KEEP_ROW) 
  table(temp_TAB$KEEP_ROW, useNA = "always")
  
  temp_KEEP_ROW_db1 <- subset(temp_TAB, temp_TAB$KEEP_ROW==1)
  temp_KEEP_ROW_db2 <- subset(temp_TAB, temp_TAB$KEEP_ROW==2)
  remove(temp_TAB)
  
  DATABASE_1 <- merge(DATABASE_1, temp_KEEP_ROW_db1[,c(1:3)], by = c("YEAR","NUTS2_CODE","CROP_CODE"), all.x = FALSE, all.y = TRUE)
  DATABASE_2 <- merge(DATABASE_2, temp_KEEP_ROW_db2[,c(1:3)], by = c("YEAR","NUTS2_CODE","CROP_CODE"), all.x = FALSE, all.y = TRUE)
  remove(temp_KEEP_ROW_db1, temp_KEEP_ROW_db2)
  
  names(DATABASE_1)==names(DATABASE_2)
  DATABASE <- rbind(DATABASE_1, DATABASE_2)
  
  DATABASE <- DATABASE[order(DATABASE$CROP_CODE, DATABASE$NUTS2_CODE, DATABASE$YEAR),] 
  remove(DATABASE_1, DATABASE_2)
  
  DATABASE$FLAG_del_row <- 0
  for (i in 2:nrow(DATABASE)){
       print(i)
       DATABASE$FLAG_del_row[i]<- ifelse(DATABASE[i,c(1:13)]==DATABASE[i-1,c(1:13)] & DATABASE$SOURCE_CODE[i]=="[apro_cpnh1]", 1, 0)}
       remove(i)
  table(DATABASE$FLAG_del_row)    #test, no repeated lines between Database_1 & Database_2

  DATABASE <- subset(DATABASE, DATABASE$NUTS2_CODE!="CY0" & DATABASE$NUTS2_CODE!="CY00" &
                               DATABASE$NUTS2_CODE!="CYZ" & DATABASE$NUTS2_CODE!="CYZZ")
  
  DATABASE <- subset(DATABASE, DATABASE$CROP_CODE!="C1100"      & DATABASE$CROP_CODE!="C1111"      &
                               DATABASE$CROP_CODE!="C1112"      & DATABASE$CROP_CODE!="C1200"      &
                               DATABASE$CROP_CODE!="C1310"      & DATABASE$CROP_CODE!="C1320"      &
                               DATABASE$CROP_CODE!="C1400"      & DATABASE$CROP_CODE!="F1111"      &
                               DATABASE$CROP_CODE!="F1112"      & DATABASE$CROP_CODE!="F1121"      &
                               DATABASE$CROP_CODE!="F1122"      & DATABASE$CROP_CODE!="F1210_1220" &
                               DATABASE$CROP_CODE!="F1212_1222" & DATABASE$CROP_CODE!="I1110-1130" &
                               DATABASE$CROP_CODE!="I1111"      & DATABASE$CROP_CODE!="I1112"      &
                               DATABASE$CROP_CODE!="S0000S"     & DATABASE$CROP_CODE!="T1100"      &
                               DATABASE$CROP_CODE!="T1200"      & DATABASE$CROP_CODE!="T1300"      &
                               DATABASE$CROP_CODE!="T1900"      & DATABASE$CROP_CODE!="T3000"      &
                               DATABASE$CROP_CODE!="V2300S"     & DATABASE$CROP_CODE!="V2710"      &
                               DATABASE$CROP_CODE!="V2720"      & DATABASE$CROP_CODE!="V3100S"     &
                               DATABASE$CROP_CODE!="V3110"      & DATABASE$CROP_CODE!="V3120"      &
                               DATABASE$CROP_CODE!="V3200S"     & DATABASE$CROP_CODE!="W1110"      &
                               DATABASE$CROP_CODE!="W1110_1120" & DATABASE$CROP_CODE!="W1120")
  
  DATABASE$rem_ROW <- grepl("Extra-Regio NUTS", DATABASE$NUTS2_NAME)
  table(DATABASE$rem_ROW)
  DATABASE <- subset(DATABASE, DATABASE$rem_ROW==FALSE)
  
#5.IMPORT ADMINISTRATIVE AREA (BY NUTS0 AND NUTS2) & CALCULATION OF AVERAGE (5 YR) ===================================================================
  TOTAL_AREA <- read.xlsx(xlsxFile = "C:\\Users\\......\\Area by NUTS 3 region - demo_r_d3area__custom_6820433_spreadsheet.xlsx", 
                          rows = c(10:1979), cols = c(1,2,11), sheet = "Sheet 1")
  colnames(TOTAL_AREA)[1:3] <- c("NUTS2_CODE", "NUTS2_NAME", "TOTAL_AREA")
  TOTAL_AREA$TOTAL_AREA <- as.numeric(TOTAL_AREA$TOTAL_AREA)
  TOTAL_AREA$TOTAL_AREA <- TOTAL_AREA$TOTAL_AREA*0.1
  
  #NUTS2 codes update - France
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR","FR",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR10","FR10",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR21","FRF2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR22","FRE2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR23","FRD2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR24","FRB0",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR25","FRD1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR26","FRC1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR30","FRE1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR41","FRF3",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR42","FRF1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR43","FRC2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR51","FRG0",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR52","FRH0",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR53","FRI3",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR61","FRI1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR62","FRJ2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR63","FRI2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR71","FRK2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR72","FRK1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR81","FRJ1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR82","FRL0",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FR83","FRM0",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FRA1","FRY1",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FRA2","FRY2",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FRA3","FRY3",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FRA4","FRY4",TOTAL_AREA$NUTS2_CODE)
  TOTAL_AREA$NUTS2_CODE <- ifelse(TOTAL_AREA$NUTS2_CODE=="FRA5","FRY5",TOTAL_AREA$NUTS2_CODE)
 
  DATABASE <- merge(DATABASE, TOTAL_AREA[,c(1,3)], by = "NUTS2_CODE", all.x = TRUE, all.y = FALSE)

  CTRL_1 <- DATABASE[is.na(DATABASE$TOTAL_AREA)==TRUE,]  
  remove(CTRL_1)

  DATABASE$LEVEL_NUTS <- nchar(DATABASE$NUTS2_CODE)
  DATABASE$LEVEL_NUTS [DATABASE$LEVEL_NUTS==2] <- "NUTS0"
  DATABASE$LEVEL_NUTS [DATABASE$LEVEL_NUTS==3] <- "NUTS1"
  DATABASE$LEVEL_NUTS [DATABASE$LEVEL_NUTS==4] <- "NUTS2"
  table(DATABASE$LEVEL_NUTS)
  
  DATABASE$LEVEL_NUTS [DATABASE$NUTS2_CODE=="EL41_42"] <- "NUTS1"
  table(DATABASE$LEVEL_NUTS)
  
  DATABASE <- subset(DATABASE, select = c("YEAR","NUTS0_CODE","NUTS0_NAME", "NUTS2_NAME","NUTS2_CODE","LEVEL_NUTS",
                                          "CROP_CODE","CROP_NAME","AREA","FLAG_AREA","PROD","FLAG_PROD","SAU","FLAG_SAU",
                                          "TOTAL_AREA","SOURCE_CODE","SOURCE_LINK","SOURCE_NAME")) 

  #Integration in the Database of 5yr average and SPA16 values 
  temp_avg <- ddply(DATABASE, c("NUTS2_CODE", "CROP_CODE"), summarise,
                    AREA_5YR_AVG = round(mean(AREA, na.rm = TRUE), 2),
                    FREQ_AREA_5YR_AVG = sum(is.na(AREA)==FALSE),
                    PROD_5YR_AVG = round(mean(PROD, na.rm = TRUE), 2),
                    FREQ_PROD_5YR_AVG = sum(is.na(PROD)==FALSE),
                    SAU_5YR_AVG  = round(mean(SAU, na.rm = TRUE), 2),
                    FREQ_SAU_5YR_AVG = sum(is.na(SAU)==FALSE))
  
  temp_avg$FREQ_AREA_5YR_AVG[is.nan(temp_avg$AREA_5YR_AVG)] <- NA
  temp_avg$AREA_5YR_AVG[is.nan(temp_avg$AREA_5YR_AVG)] <- NA
  
  temp_avg$FREQ_PROD_5YR_AVG[is.nan(temp_avg$PROD_5YR_AVG)] <- NA
  temp_avg$PROD_5YR_AVG[is.nan(temp_avg$PROD_5YR_AVG)] <- NA
  
  temp_avg$FREQ_SAU_5YR_AVG[is.nan(temp_avg$SAU_5YR_AVG)] <- NA
  temp_avg$SAU_5YR_AVG[is.nan(temp_avg$SAU_5YR_AVG)] <- NA
  
  DATABASE <- merge(DATABASE, temp_avg, by = c("NUTS2_CODE", "CROP_CODE"), all.x = TRUE, all.y = FALSE)  
  remove(temp_avg)
  
  DATABASE$PERC_AREA_SAU <- round(DATABASE$AREA/DATABASE$SAU*100, 2)
  DATABASE$PERC_AREA_totAREA <- round(DATABASE$AREA/DATABASE$TOTAL_AREA*100, 2)
 
  DATABASE <- subset(DATABASE, select = c("YEAR","NUTS0_CODE","NUTS0_NAME", "NUTS2_NAME","NUTS2_CODE","LEVEL_NUTS",
                                          "CROP_CODE","CROP_NAME","AREA","FLAG_AREA","PROD","FLAG_PROD","SAU","FLAG_SAU",
                                          "TOTAL_AREA","AREA_5YR_AVG","FREQ_AREA_5YR_AVG","PROD_5YR_AVG","FREQ_PROD_5YR_AVG",
                                          "SAU_5YR_AVG","FREQ_SAU_5YR_AVG","PERC_AREA_SAU","PERC_AREA_totAREA",
                                          "SOURCE_CODE","SOURCE_LINK","SOURCE_NAME")) 
  
#6.CALCULATION SPA16 ESTIMATES ====================================================================================================================
  options(scipen = 999)
  
  setwd("C:\\Users\\................\\data folder")
  DB_EU_SPA16 <-  read.xlsx(xlsxFile = "Dati SPA 2016 - EU27 (1502_23).xlsx", rows = c(12:342), cols = c(1:75), sheet = "Sheet 1")
  
  for (i in 3:ncol(DB_EU_SPA16)) {
       DB_EU_SPA16[,i] <- as.numeric(DB_EU_SPA16[,i])}
       remove(i)
  
  DB_EU_SPA16_PERC <- subset(DB_EU_SPA16, select = c("NUTS2_CODE","NUTS2_NAME",
                                                     "F1000","F2000","F3000","F4000",
                                                     "W1110_1120","W1190","W1200","W1300",
                                                     "O1100","O1910"))
  DB_EU_SPA16_PERC$NUTS0_CODE <- substr(DB_EU_SPA16_PERC$NUTS2_CODE, 1, 2)
  
  temp_TAB_NUTS0 <- data.frame(NUTS0_CODE = DB_EU_SPA16$NUTS2_CODE, NUTS0_NAME = DB_EU_SPA16$NUTS2_NAME)
  temp_TAB_NUTS0 <- subset(temp_TAB_NUTS0, nchar(temp_TAB_NUTS0$NUTS0_CODE)==2)
  DB_EU_SPA16_PERC <- merge(DB_EU_SPA16_PERC, temp_TAB_NUTS0, by = c("NUTS0_CODE"), all.x = TRUE, all.y = FALSE)
  remove(temp_TAB_NUTS0, DB_EU_SPA16)  
  
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[,c(1,14,2:13)]
  
  #NUTS2 codes update - France
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR","FR",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR10","FR10",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR21","FRF2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR22","FRE2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR23","FRD2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR24","FRB0",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR25","FRD1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR26","FRC1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR30","FRE1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR41","FRF3",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR42","FRF1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR43","FRC2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR51","FRG0",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR52","FRH0",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR53","FRI3",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR61","FRI1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR62","FRJ2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR63","FRI2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR71","FRK2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR72","FRK1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR81","FRJ1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR82","FRL0",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FR83","FRM0",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FRA1","FRY1",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FRA2","FRY2",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FRA3","FRY3",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_CODE <- ifelse(DB_EU_SPA16_PERC$NUTS2_CODE=="FRA4","FRY4",DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="France","France",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Île de France","Île de France",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Champagne-Ardenne (NUTS 2013)","Champagne-Ardenne",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Picardie (NUTS 2013)","Picardie",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Haute-Normandie (NUTS 2013)","Haute-Normandie",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Centre (FR) (NUTS 2013)","Centre - Val de Loire",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Basse-Normandie (NUTS 2013)","Basse-Normandie",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Bourgogne (NUTS 2013)","Bourgogne",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Nord-Pas-de-Calais (NUTS 2013)","Nord-Pas-de-Calais",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Lorraine (NUTS 2013)","Lorraine",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Alsace (NUTS 2013)","Alsace",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Franche-Comté (NUTS 2013)","Franche-Comté",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Pays de la Loire (NUTS 2013)","Pays-de-la-Loire",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Bretagne (NUTS 2013)","Bretagne",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Poitou-Charentes (NUTS 2013)","Poitou-Charentes",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Aquitaine (NUTS 2013)","Aquitaine",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Midi-Pyrénées (NUTS 2013)","Midi-Pyrénées",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Limousin (NUTS 2013)","Limousin",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Rhône-Alpes (NUTS 2013)","Rhône-Alpes",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Auvergne (NUTS 2013)","Auvergne",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Languedoc-Roussillon (NUTS 2013)","Languedoc-Roussillon",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Provence-Alpes-Côte d'Azur (NUTS 2013)","Provence-Alpes-Côte d'Azur",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Corse (NUTS 2013)","Corse",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Guadeloupe (NUTS 2013)","Guadeloupe",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Martinique (NUTS 2013)","Martinique",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="Guyane (NUTS 2013)","Guyane",DB_EU_SPA16_PERC$NUTS2_NAME)
  DB_EU_SPA16_PERC$NUTS2_NAME <- ifelse(DB_EU_SPA16_PERC$NUTS2_NAME=="La Réunion (NUTS 2013)","La Réunion",DB_EU_SPA16_PERC$NUTS2_NAME)
  
  ADD_Mayotte_FR <- data.frame("FR","France","FRY5","Mayotte",0,0,0,0,0,0,0,0,0,0)
  colnames(ADD_Mayotte_FR)[1:14] <- colnames(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC <- rbind(DB_EU_SPA16_PERC, ADD_Mayotte_FR)
  remove(ADD_Mayotte_FR)
  
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[order(DB_EU_SPA16_PERC$NUTS2_CODE),]
  DB_EU_SPA16_PERC[is.na(DB_EU_SPA16_PERC)] <- 0 
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC <- subset(DB_EU_SPA16_PERC, NUTS2_CODE!="FR91" & NUTS2_CODE!="FR92"&NUTS2_CODE!="FR93" & NUTS2_CODE!="FR94")
  detach(DB_EU_SPA16_PERC)
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC$F0000 <- F1000+F2000+F3000+F4000
  DB_EU_SPA16_PERC$W1000 <- W1110_1120+W1200+W1300+W1190
  DB_EU_SPA16_PERC$O1000 <- O1100+O1910 
  detach(DB_EU_SPA16_PERC)
  
  #Calculation of the incidence of crop area (AREA) on national Utilised agricultural area (UAA)
  TEMP_sup_naz <- subset(DB_EU_SPA16_PERC[,c(1:2,15:17)], nchar(DB_EU_SPA16_PERC$NUTS2_CODE)==2)
  colnames(TEMP_sup_naz)[3:5] <- c("F0000_naz","W1000_naz","O1000_naz")
  DB_EU_SPA16_PERC <- merge(DB_EU_SPA16_PERC, TEMP_sup_naz[,c(1,3:5)], by = "NUTS0_CODE", all.x = TRUE, all.y = FALSE)
  remove(TEMP_sup_naz)
  
  CODE_EU_27 <- c("BE",	"BG",	"CZ",	"DK",	"DE",	"EE",	"IE",	"EL",	"ES",	"FR",	"HR",	"IT",	
                  "CY",	"LV",	"LT",	"LU",	"HU",	"MT",	"NL",	"AT",	"PL",	"PT",	"RO",	"SI",	
                  "SK",	"FI",	"SE")
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[DB_EU_SPA16_PERC$NUTS0_CODE %in% CODE_EU_27,]
  remove(CODE_EU_27) 
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC$PERC_F0000 <- F0000/F0000_naz
  DB_EU_SPA16_PERC$PERC_F1000 <- F1000/F0000_naz
  DB_EU_SPA16_PERC$PERC_F2000 <- F2000/F0000_naz
  DB_EU_SPA16_PERC$PERC_F3000 <- F3000/F0000_naz
  DB_EU_SPA16_PERC$PERC_F4000 <- F4000/F0000_naz
  
  DB_EU_SPA16_PERC$PERC_W1000      <- W1000/W1000_naz
  DB_EU_SPA16_PERC$PERC_W1110_1120 <- W1110_1120/W1000_naz
  DB_EU_SPA16_PERC$PERC_W1190      <- W1190/W1000_naz
  DB_EU_SPA16_PERC$PERC_W1200      <- W1200/W1000_naz
  DB_EU_SPA16_PERC$PERC_W1300      <- W1300/W1000_naz
  
  DB_EU_SPA16_PERC$PERC_O1000 <- O1000/O1000_naz
  DB_EU_SPA16_PERC$PERC_O1100 <- O1100/O1000_naz
  DB_EU_SPA16_PERC$PERC_O1910 <- O1910/O1000_naz
  detach(DB_EU_SPA16_PERC)
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC$PERC_F0000      <- ifelse(F0000_naz==0 & is.na(PERC_F0000)==TRUE     , 0, PERC_F0000)
  DB_EU_SPA16_PERC$PERC_F1000      <- ifelse(F0000_naz==0 & is.na(PERC_F1000)==TRUE     , 0, PERC_F1000)
  DB_EU_SPA16_PERC$PERC_F2000      <- ifelse(F0000_naz==0 & is.na(PERC_F2000)==TRUE     , 0, PERC_F2000)
  DB_EU_SPA16_PERC$PERC_F3000      <- ifelse(F0000_naz==0 & is.na(PERC_F3000)==TRUE     , 0, PERC_F3000)
  DB_EU_SPA16_PERC$PERC_F4000      <- ifelse(F0000_naz==0 & is.na(PERC_F4000)==TRUE     , 0, PERC_F4000)
  DB_EU_SPA16_PERC$PERC_W1000      <- ifelse(W1000_naz==0 & is.na(PERC_W1000)==TRUE     , 0, PERC_W1000)
  DB_EU_SPA16_PERC$PERC_W1110_1120 <- ifelse(W1000_naz==0 & is.na(PERC_W1110_1120)==TRUE, 0, PERC_W1110_1120)
  DB_EU_SPA16_PERC$PERC_W1190      <- ifelse(W1000_naz==0 & is.na(PERC_W1190)==TRUE     , 0, PERC_W1190)
  DB_EU_SPA16_PERC$PERC_W1200      <- ifelse(W1000_naz==0 & is.na(PERC_W1200)==TRUE     , 0, PERC_W1200)
  DB_EU_SPA16_PERC$PERC_W1300      <- ifelse(W1000_naz==0 & is.na(PERC_W1300)==TRUE     , 0, PERC_W1300)
  DB_EU_SPA16_PERC$PERC_O1000      <- ifelse(O1000_naz==0 & is.na(PERC_O1000)==TRUE     , 0, PERC_O1000)
  DB_EU_SPA16_PERC$PERC_O1100      <- ifelse(O1000_naz==0 & is.na(PERC_O1100)==TRUE     , 0, PERC_O1100)
  DB_EU_SPA16_PERC$PERC_O1910      <- ifelse(O1000_naz==0 & is.na(PERC_O1910)==TRUE     , 0, PERC_O1910)
  detach(DB_EU_SPA16_PERC)
  
  DB_EU_SPA16_PERC[sapply(DB_EU_SPA16_PERC, is.infinite)] <- NA
  DB_EU_SPA16_PERC[sapply(DB_EU_SPA16_PERC, is.nan)] <- NA
  
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[,c(1:4,21:33)]
  DB_EU_SPA16_PERC <- melt(DB_EU_SPA16_PERC, id = c("NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME"))
  colnames(DB_EU_SPA16_PERC)[5:6] <- c("CROP_CODE","PERC_SPA16")
  DB_EU_SPA16_PERC$CROP_CODE <- gsub("PERC_", "", DB_EU_SPA16_PERC$CROP_CODE)
  
  DB_EU_SPA16_PERC$CROP_CODE_2[substr(DB_EU_SPA16_PERC$CROP_CODE,1,1)=="F"] <- "F0000"
  DB_EU_SPA16_PERC$CROP_CODE_2[substr(DB_EU_SPA16_PERC$CROP_CODE,1,1)=="O"] <- "O1000"
  DB_EU_SPA16_PERC$CROP_CODE_2[substr(DB_EU_SPA16_PERC$CROP_CODE,1,1)=="W"] <- "W1000"
  table(DB_EU_SPA16_PERC$CROP_CODE_2)
  
  TEMP_17 <- data.frame(YEAR = 2017, DB_EU_SPA16_PERC)
  TEMP_18 <- data.frame(YEAR = 2018, DB_EU_SPA16_PERC)
  TEMP_19 <- data.frame(YEAR = 2019, DB_EU_SPA16_PERC)
  TEMP_20 <- data.frame(YEAR = 2020, DB_EU_SPA16_PERC)
  TEMP_21 <- data.frame(YEAR = 2021, DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC <- rbind(TEMP_17, TEMP_18, TEMP_19, TEMP_20, TEMP_21)
  remove(TEMP_17, TEMP_18, TEMP_19, TEMP_20, TEMP_21)
  
  CODE_EU_27 <- c("BE",	"BG",	"CZ",	"DK",	"DE",	"EE",	"IE",	"EL",	"ES",	"FR",	"HR",	"IT",	
                  "CY",	"LV",	"LT",	"LU",	"HU",	"MT",	"NL",	"AT",	"PL",	"PT",	"RO",	"SI",	
                  "SK",	"FI",	"SE")
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[DB_EU_SPA16_PERC$NUTS0_CODE %in% CODE_EU_27,]
  colnames(DB_EU_SPA16_PERC)[8] <- "CROP_CODE_MACRO_COLT"
  remove(CODE_EU_27)
  
  temp_DB <- subset(DATABASE, select = c("YEAR","NUTS0_CODE","NUTS2_CODE","CROP_CODE","AREA",
                                         "PROD","SAU"))
  
  temp_sup_naz <- subset(temp_DB, temp_DB$CROP_CODE=="F0000" | temp_DB$CROP_CODE=="O1000" | temp_DB$CROP_CODE=="W1000")
  temp_sup_naz <- subset(temp_sup_naz, nchar(temp_sup_naz$NUTS2_CODE)==2)
  temp_sup_naz <- subset(temp_sup_naz, is.na(temp_sup_naz$AREA)==FALSE & is.na(temp_sup_naz$PROD)==FALSE)
  
  DB_EU_SPA16_PERC <- merge(DB_EU_SPA16_PERC, temp_sup_naz[,c(1,2,4,5:7)], by.x = c("YEAR","NUTS0_CODE","CROP_CODE_MACRO_COLT"), by.y = c("YEAR","NUTS0_CODE","CROP_CODE"), all.x = TRUE, all.y = FALSE)
  
  #Crop estimate calculation SPA16
  DB_EU_SPA16_PERC$AREA_estmt_SPA16 <- DB_EU_SPA16_PERC$AREA*DB_EU_SPA16_PERC$PERC_SPA16
  DB_EU_SPA16_PERC$PROD_estmt_SPA16 <- DB_EU_SPA16_PERC$PROD*DB_EU_SPA16_PERC$PERC_SPA16
  DB_EU_SPA16_PERC$SAU_estmt_SPA16 <- DB_EU_SPA16_PERC$SAU*DB_EU_SPA16_PERC$PERC_SPA16
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[,c(1:7,12:14)]
  colnames(DB_EU_SPA16_PERC)[8:10] <- c("AREA","PROD","SAU")
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC$FLAG_AREA <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "estmt_FSS16", "")
  DB_EU_SPA16_PERC$FLAG_PROD <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "estmt_FSS16", "")
  DB_EU_SPA16_PERC$FLAG_SAU  <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "estmt_FSS16", "")
  DB_EU_SPA16_PERC <- DB_EU_SPA16_PERC[order(DB_EU_SPA16_PERC$NUTS0_CODE, DB_EU_SPA16_PERC$YEAR, DB_EU_SPA16_PERC$CROP_CODE),]
  detach(DB_EU_SPA16_PERC)
  
  DB_EU_SPA16_PERC$CROP_NAME <- ""
  
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="F0000",	"Fruits, berries and nuts (excluding citrus fruits, grapes and strawberries [F0000])", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="W1000",	"Grapes [W1000]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="O1000",	"Olives [O1000]", DB_EU_SPA16_PERC$CROP_NAME)
  
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="F1000",	"Fruits from temperate climate zones [F1000]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="F2000",	"Fruits from subtropical and tropical climate zones [F2000]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="F3000",	"Berries (excluding strawberries) [F3000]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="F4000",	"Nuts [F4000]", DB_EU_SPA16_PERC$CROP_NAME)
  
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="W1110_1120",	"Grapes for wines with geographical indication (PDO/PGI) [W1110_1120]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="W1190",	"Grapes for other wines n.e.c. (without PDO/PGI) [W1190]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="W1200",	"Grapes for table use [W1200]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="W1300",	"Grapes for raisins [W1300]", DB_EU_SPA16_PERC$CROP_NAME)
  
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="O1100",	"Olives for table use [O1100]", DB_EU_SPA16_PERC$CROP_NAME)
  DB_EU_SPA16_PERC$CROP_NAME <- ifelse(DB_EU_SPA16_PERC$CROP_CODE=="O1910",	"Olives for oil [O1910]", DB_EU_SPA16_PERC$CROP_NAME)
  
  table(DB_EU_SPA16_PERC$CROP_NAME)  
  
  DB_EU_SPA16_PERC$CROP_CODE_MACRO <- NULL
  remove(temp_sup_naz)
  
  #Calculation of average 5yr SPA16 data and Database integration 
  temp_avg <- ddply(DB_EU_SPA16_PERC, c("NUTS2_CODE", "CROP_CODE"), summarise,
                    AREA_5YR_AVG = round(mean(AREA, na.rm = TRUE), 2),
                    FREQ_AREA_5YR_AVG = sum(is.na(AREA)==FALSE),
                    PROD_5YR_AVG = round(mean(PROD, na.rm = TRUE), 2),
                    FREQ_PROD_5YR_AVG = sum(is.na(PROD)==FALSE),
                    SAU_5YR_AVG  = round(mean(SAU, na.rm = TRUE), 2),
                    FREQ_SAU_5YR_AVG = sum(is.na(SAU)==FALSE))
  
  temp_avg$FREQ_AREA_5YR_AVG[is.nan(temp_avg$AREA_5YR_AVG)] <- NA
  temp_avg$AREA_5YR_AVG[is.nan(temp_avg$AREA_5YR_AVG)] <- NA
  
  temp_avg$FREQ_PROD_5YR_AVG[is.nan(temp_avg$PROD_5YR_AVG)] <- NA
  temp_avg$PROD_5YR_AVG[is.nan(temp_avg$PROD_5YR_AVG)] <- NA
  
  temp_avg$FREQ_SAU_5YR_AVG[is.nan(temp_avg$SAU_5YR_AVG)] <- NA
  temp_avg$SAU_5YR_AVG[is.nan(temp_avg$SAU_5YR_AVG)] <- NA
  
  DB_EU_SPA16_PERC <- merge(DB_EU_SPA16_PERC, temp_avg, by = c("NUTS2_CODE", "CROP_CODE"), all.x = TRUE, all.y = FALSE)  
  remove(temp_avg)
  
  DB_EU_SPA16_PERC$PERC_AREA_SAU <- round(DB_EU_SPA16_PERC$AREA/DB_EU_SPA16_PERC$SAU*100, 2)
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC$SOURCE_CODE <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "[EF_LUS_ALLCROPS]","")
  DB_EU_SPA16_PERC$SOURCE_LINK <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "https://ec.europa.eu/eurostat/databrowser/view/EF_LUS_ALLCROPS__custom_4951867/default/table?lang=en","")
  DB_EU_SPA16_PERC$SOURCE_NAME <- ifelse(CROP_CODE!="F0000" & CROP_CODE!="O1000" & CROP_CODE!="W1000", "Crops by classes of utilised agricultural area in number of farms and hectare by NUTS 2 regions","")
  detach(DB_EU_SPA16_PERC)
  
  attach(DB_EU_SPA16_PERC)
  DB_EU_SPA16_PERC <- subset(DB_EU_SPA16_PERC, SOURCE_NAME!="" & SOURCE_LINK!="" & SOURCE_CODE!="")
  detach(DB_EU_SPA16_PERC)
  
  #Matching Administrative AREA to SPA16 
  DB_EU_SPA16_PERC <- merge(DB_EU_SPA16_PERC, TOTAL_AREA[,c(1,3)], by = "NUTS2_CODE", all.x = TRUE, all.y = FALSE)
  
  DB_EU_SPA16_PERC$LEVEL_NUTS <- nchar(DB_EU_SPA16_PERC$NUTS2_CODE)
  DB_EU_SPA16_PERC$LEVEL_NUTS [DB_EU_SPA16_PERC$LEVEL_NUTS==2] <- "NUTS0"
  DB_EU_SPA16_PERC$LEVEL_NUTS [DB_EU_SPA16_PERC$LEVEL_NUTS==3] <- "NUTS1"
  DB_EU_SPA16_PERC$LEVEL_NUTS [DB_EU_SPA16_PERC$LEVEL_NUTS==4] <- "NUTS2"
  table(DB_EU_SPA16_PERC$LEVEL_NUTS)
  
  colnames(DB_EU_SPA16_PERC)
  colnames(DATABASE)
  DB_EU_SPA16_PERC$PERC_AREA_totAREA <- round((DB_EU_SPA16_PERC$AREA/DB_EU_SPA16_PERC$TOTAL_AREA)*100, 5)
  remove(TOTAL_AREA)
  
#7.SPA_2016 ALIGNMENT FOR APPEND IN DATABASE =======================================================================================================
  DB_EU_SPA16_PERC <- subset(DB_EU_SPA16_PERC, select = colnames(DATABASE))  
  
  names(DATABASE) == names(DB_EU_SPA16_PERC)  
  
  DATABASE_final <- rbind(DATABASE, DB_EU_SPA16_PERC)
  
  DATABASE_final <- DATABASE_final[order(DATABASE_final$NUTS2_CODE, DATABASE_final$CROP_CODE, DATABASE_final$YEAR),]
  
  DATABASE_final <- merge(DATABASE_final, TEMP[,c(1,3:7)], by = c("NUTS2_CODE","CROP_CODE"), all.x = TRUE, all.y = FALSE)
  
  DATABASE_final <- subset(DATABASE_final, select = c("YEAR","NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME","LEVEL_NUTS","CROP_CODE","CROP_NAME",       
                                                      "AREA","FLAG_AREA","AREA_5YR_AVG","FREQ_AREA_5YR_AVG",
                                                      "PROD","FLAG_PROD","PROD_5YR_AVG","FREQ_PROD_5YR_AVG",    
                                                      "SAU","FLAG_SAU","SAU_5YR_AVG","FREQ_SAU_5YR_AVG",
                                                      "PERC_AREA_SAU","MEAN_PERC_AREA_SAU","FREQ_PERC_AREA_SAU",
                                                      "TOTAL_AREA","PERC_AREA_totAREA","MEAN_PERC_AREA_totAREA","FREQ_PERC_AREA_totAREA",
                                                      "SOURCE_CODE","SOURCE_LINK","SOURCE_NAME"))
  colnames(DATABASE_final)[1:30] <- c("YEAR","NUTS0_CODE","NUTS0_NAME","NUTS2_CODE","NUTS2_NAME","NUTS_LEVEL","CROP_CODE","CROP_NAME",
                                      "AREA","AREA_FLG","AREA_AV","AREA_FRQ_AV",
                                      "PROD","PROD_FLG","PROD_AV","PROD_FRQ_AV",
                                      "UAA","UAA_FLG","UAA_AV","UAA_FRQ_AV",
                                      "PCR_UAA","PCR_UAA_AV","PCR_UAA_FRQ_AV",
                                      "NUTS_AREA","PCR_NUTS_AREA","PCR_NUTS_AV","PCR_NUTS_FRQ_AV",
                                      "SOURCE_CODE","SOURCE_LINK","SOURCE_NAME")
                
  #Eliminate duplicate rows between spa16 and Database
  DATABASE_final <- DATABASE_final[order(DATABASE_final$NUTS2_CODE, DATABASE_final$CROP_CODE, DATABASE_final$YEAR),]
  
  DATABASE_final$FLAG_del_row <- 0
  for (i in 2:nrow(DATABASE_final)){
       print(i)
       DATABASE_final$FLAG_del_row[i]<- ifelse(DATABASE_final[i,c(1:8)]==DATABASE_final[i-1,c(1:8)] & DATABASE_final$SOURCE_CODE[i]=="[EF_LUS_ALLCROPS]", 1, 0)}
       remove(i)
       
  table(DATABASE_final$FLAG_del_row)
  DATABASE_final <- subset(DATABASE_final, DATABASE_final$FLAG_del_row!=1)
  
  DATABASE_final$PCR_UAA_AV[is.nan(DATABASE_final$PCR_UAA_AV)] <- NA
  DATABASE_final$PCR_NUTS_AV[is.nan(DATABASE_final$PCR_NUTS_AV)] <- NA
  
  TAB_FREQ_DB_final <- ddply(DATABASE_final, c("NUTS2_CODE","CROP_CODE"), summarise,
                             FREQ_AREA = sum(is.na(AREA)==FALSE),
                             FREQ_PROD = sum(is.na(PROD)==FALSE),
                             FREQ_UAA  = sum(is.na(UAA)==FALSE))
  
  table(DATABASE_final$NUTS0_CODE)
  
  TAB_FREQ_NUTS2_CODE <- ddply(unique(DATABASE_final[c(2,4)]), c("NUTS0_CODE"), summarise,
                               FREQ_NUTS2  = length(is.na(NUTS2_CODE)==FALSE))
  
  colnames(DATABASE_final)
  DATABASE_final$FLAG_del_row <- NULL
  
  remove(temp_DB)
  
  #Export final Database Crops
  setwd("C:\\Users\\......\\Data_folder")
  wb <- createWorkbook()
  addWorksheet(wb, sheetName = "DATABASE_Final_Release")
  writeData(wb, sheet = 1, x = DATABASE_final, na.string = NA)
  saveWorkbook(wb, "Database CROPS Eurostat - Final Release.xlsx", overwrite = TRUE)
  remove(wb)
  
  remove(temp_input_CTRL, temp_CTRL_freq_all_nat)