rm(list=ls())
setwd("E:/Dropbox/Research/_DeLoach/Ag_Data_News/")

#############################################################
## Code to create plots in Ag Data News article            ##
##  Title: What is a Farm?       ##
##  Link:  https://asmith.ucdavis.edu/news//what-farm    ##
#############################################################

# If a package is installed, it will be loaded. If any are not, the missing package(s) will be installed 
# from CRAN and then loaded.
pacman::p_load(Quandl,tidyverse,readxl,RColorBrewer,httr,blscrapeR,rnassqs,data.table,tigris,tmap,pals,ggpubr)


########################################################
##  Farms
########################################################

# api for NASS (get yours at quickstats.nass.usda.gov/api)
apikey <- "INSERT YOUR API KEY HERE"

nassqs_auth(apikey) # authentificate your api key
nassqs_params() # check what variables are there
domain_list <- nassqs_param_values("domain_desc") # get commodity list


# set years
# you will receive an "bad request - invalid query" error if there are no data for your chosen years.
years <- c(2002,2007,2012,2017)

# set parameters for API calls except for years 
# (we are not allowed to request many records at once, so I iterate API calls by years)
params <-list(source_desc="CENSUS",domain_desc="ECONOMIC CLASS",agg_level_desc="NATIONAL") 

# create data based on the parameters and years using 'nassqs'
df <- map(years, function(x){params[["year"]] <- x; nassqs(params)}) %>% bind_rows() %>% as.data.table()


##############################################################
# Plot Number of Operations by Economic Class
##############################################################
plot_bar_operations <- df %>%
  filter(short_desc == "COMMODITY TOTALS - OPERATIONS WITH SALES") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (1,000,000 TO 2,499,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (2,500,000 TO 4,999,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (5,000,000 OR MORE $)") %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (LESS THAN 1,000 $)","<1",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000 TO 2,499 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (2,500 TO 4,999 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (5,000 TO 9,999 $)","5-10",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (10,000 TO 24,999 $)","10-25",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (25,000 TO 49,999 $)","25-50",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (50,000 TO 99,999 $)","50-100",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (100,000 TO 249,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (250,000 TO 499,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (500,000 TO 999,999 $)","500-1,000",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000,000 OR MORE $)",">1,000",domaincat_desc)) %>%
  mutate(Value=as.numeric(gsub(",","",Value))/1000) %>%
  select(domaincat_desc,year,Value)%>%
  group_by(domaincat_desc,year) %>%
  summarise(Value=sum(Value)) %>%
  ggplot(aes(fill=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
)),y=Value,x=factor(year,levels=c(2017,2012,2007,2002)),label=floor(Value))) +
    geom_bar(position="stack",stat="identity",width=.9) +
    geom_text(size = 3, position = position_stack(vjust = .5))+
    ggtitle(paste0("Thousands of Farms at Various Sales Levels"))+
    coord_flip() +
    theme_minimal()+
    scale_fill_brewer(palette="Spectral")+
    labs(x=NULL,y="Thousands of Farms",fill="Actual Sales ($000)")+
    theme(plot.title = element_text(hjust = 0.5,size=16),legend.text = element_text(hjust = 1,size=12), text = element_text(size=14))

# Draw and save plot
plot_bar_operations
ggsave("farm_operations.png")



##############################################################
# Plot Number of Operations, Sales and Land by Economic Class
##############################################################
plot_bar_operations2 <- df %>%
  filter(short_desc %in% c("COMMODITY TOTALS - OPERATIONS WITH SALES","COMMODITY TOTALS - SALES, MEASURED IN $","AG LAND, CROPLAND, HARVESTED - ACRES","AG LAND, CROPLAND, (EXCL HARVESTED & PASTURED) - ACRES","AG LAND, PASTURELAND - ACRES","AG LAND, WOODLAND, (EXCL PASTURED) - ACRES","AG LAND, (EXCL CROPLAND & PASTURELAND & WOODLAND) - ACRES")) %>%
  mutate(short_desc=ifelse(short_desc %in% c("AG LAND, CROPLAND, HARVESTED - ACRES","AG LAND, CROPLAND, (EXCL HARVESTED & PASTURED) - ACRES","AG LAND, PASTURELAND - ACRES","AG LAND, WOODLAND, (EXCL PASTURED) - ACRES","AG LAND, (EXCL CROPLAND & PASTURELAND & WOODLAND) - ACRES"),"Land",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="COMMODITY TOTALS - OPERATIONS WITH SALES","Farms",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="COMMODITY TOTALS - SALES, MEASURED IN $","Sales",short_desc)) %>%
  filter(year==2017) %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (1,000,000 TO 2,499,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (2,500,000 TO 4,999,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (5,000,000 OR MORE $)") %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (LESS THAN 1,000 $)","<1",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000 TO 2,499 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (2,500 TO 4,999 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (5,000 TO 9,999 $)","5-10",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (10,000 TO 24,999 $)","10-25",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (25,000 TO 49,999 $)","25-50",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (50,000 TO 99,999 $)","50-100",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (100,000 TO 249,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (250,000 TO 499,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (500,000 TO 999,999 $)","500-1,000",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000,000 OR MORE $)",">1,000",domaincat_desc)) %>%
  mutate(Value=as.numeric(gsub(",","",Value))/1000) %>%
  select(domaincat_desc,short_desc,Value)%>%
  group_by(domaincat_desc,short_desc) %>%
  summarise(Value=sum(Value)) %>%
  ggplot(aes(fill=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
  )),y=Value,x=factor(short_desc, levels=c("Sales","Land","Farms")))) +
    geom_bar(stat="identity",position = position_fill(),width=.9) +
    #geom_text(size = 3, position = position_stack(vjust = .5))+
    ggtitle(paste0("Farms, Land Use, and Sales in 2017"))+
    coord_flip() +
    theme_minimal()+
    scale_fill_brewer(palette="Spectral")+
    labs(x=NULL,y="Percent",fill="Actual Sales ($000)")+
    scale_y_continuous(labels = scales::percent_format())+
    theme(plot.title = element_text(hjust = 0.5,size=16),legend.text = element_text(hjust = 1,size=12), text = element_text(size=14))

# Draw and save plot
plot_bar_operations2
ggsave("farm_operations2.png")





##############################################################
# Plot Net Income by Economic Class and Whether Made or Lost Money
##############################################################
plot_bar_netincome_pl <- df %>%
  filter(short_desc %in% c("INCOME, NET CASH FARM, OF PRODUCERS - OPERATIONS WITH GAIN","INCOME, NET CASH FARM, OF PRODUCERS - OPERATIONS WITH LOSS")) %>%
  filter(year==2017) %>%
  mutate(short_desc=ifelse(short_desc=="INCOME, NET CASH FARM, OF PRODUCERS - OPERATIONS WITH GAIN","Profit","Loss")) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (LESS THAN 1,000 $)","<1",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000 TO 2,499 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (2,500 TO 4,999 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (5,000 TO 9,999 $)","5-10",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (10,000 TO 24,999 $)","10-25",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (25,000 TO 49,999 $)","25-50",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (50,000 TO 99,999 $)","50-100",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (100,000 TO 249,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (250,000 TO 499,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (500,000 TO 999,999 $)","500-1,000",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000,000 OR MORE $)",">1,000",domaincat_desc))%>%
  mutate(Value=as.numeric(gsub(",","",Value))/1000) %>%
  select(domaincat_desc,short_desc,year,Value) %>%
  group_by(domaincat_desc,short_desc) %>%
  summarise(Value=sum(Value)) %>%
  mutate(Valuelabel=ifelse(Value>20,Value,NA))%>%
  ggplot(aes(fill=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
  )),y=Value,x=factor(short_desc,levels=c("Loss","Profit")),label=floor(Valuelabel))) +
    geom_bar(position="stack",stat="identity",width=.9) +
    geom_text(size = 3, position = position_stack(vjust = .5))+
    ggtitle(paste0("Thousands of Farms by Whether they Made or Lost Money in 2017"))+
    coord_flip() +
    theme_minimal()+
    scale_fill_brewer(palette="Spectral")+
    labs(x=NULL,y="Thousands of Farms",fill="Actual Sales ($000)")+
    theme(plot.title = element_text(hjust = 0.5,size=16),legend.text = element_text(hjust = 1,size=12), text = element_text(size=14))

# Draw and save plot
plot_bar_netincome_pl
ggsave("farm_netincome_pl.png")





##############################################################
# Plot Net Income by Economic Class
##############################################################
plot_bar_netincome <- df %>%
  filter(short_desc == "INCOME, NET CASH FARM, OF OPERATIONS - NET INCOME, MEASURED IN $ / OPERATION") %>%
  filter(year==2017) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (LESS THAN 1,000 $)","<1",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000 TO 2,499 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (2,500 TO 4,999 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (5,000 TO 9,999 $)","5-10",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (10,000 TO 24,999 $)","10-25",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (25,000 TO 49,999 $)","25-50",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (50,000 TO 99,999 $)","50-100",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (100,000 TO 249,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (250,000 TO 499,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (500,000 TO 999,999 $)","500-1,000",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000,000 OR MORE $)",">1,000",domaincat_desc))%>%
  mutate(Value=as.numeric(gsub(",","",Value))/1000) %>%
  select(domaincat_desc,year,Value) %>%
  group_by(domaincat_desc,year) %>%
  summarise(Value=sum(Value)) %>%
  ggplot(aes(x=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
  )),y=Value,fill=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
  )),label=floor(Value))) +
    geom_bar(position="stack",stat="identity",width=.9,show.legend = FALSE) +
    geom_text(size = 4, position = position_stack(vjust = 1))+
    ggtitle(paste0("Average Net Farm Income per Farm by Sales Level (2017)"))+
    #coord_flip() +
    theme_minimal()+
    scale_fill_brewer(palette="Spectral")+
    labs(x="Actual Sales ($000)",y="Thousands of Dollars")+
    theme(plot.title = element_text(hjust = 0.5,size=16),text = element_text(size=14))

# Draw and save plot
plot_bar_netincome
ggsave("farm_netincome.png")




# Plot Number of Operations by Economic Class
plot_bar_activity <- df %>%
  filter(year==2017) %>%
  filter(short_desc %in% c("COMMODITY TOTALS - OPERATIONS WITH SALES","TRACTORS - OPERATIONS WITH INVENTORY","TRUCKS, INCL PICKUPS - OPERATIONS WITH INVENTORY","AG LAND, OWNED - NUMBER OF OPERATIONS","AG LAND, RENTED FROM OTHERS - NUMBER OF OPERATIONS","AG LAND, RENTED TO OTHERS - NUMBER OF OPERATIONS","FARM OPERATIONS, ORGANIZATION, TAX PURPOSES, FAMILY & INDIVIDUAL - NUMBER OF OPERATIONS","GRAIN - OPERATIONS WITH SALES","HAY & HAYLAGE - OPERATIONS WITH AREA HARVESTED","FRUIT & TREE NUT TOTALS - OPERATIONS WITH SALES","VEGETABLE TOTALS, INCL SEEDS & TRANSPLANTS, IN THE OPEN - OPERATIONS WITH SALES","CATTLE, INCL CALVES - OPERATIONS WITH INVENTORY","CHICKENS, LAYERS - OPERATIONS WITH INVENTORY","CHICKENS, PULLETS, REPLACEMENT - OPERATIONS WITH INVENTORY")
) %>%
  mutate(short_desc=ifelse(short_desc=="COMMODITY TOTALS - OPERATIONS WITH SALES","Total",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="TRACTORS - OPERATIONS WITH INVENTORY","Own Tractor",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="TRUCKS, INCL PICKUPS - OPERATIONS WITH INVENTORY","Own Truck",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="AG LAND, OWNED - NUMBER OF OPERATIONS","Own Land",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="AG LAND, RENTED FROM OTHERS - NUMBER OF OPERATIONS","Rent Land from Others",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="AG LAND, RENTED TO OTHERS - NUMBER OF OPERATIONS","Rent Land to Others",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="FARM OPERATIONS, ORGANIZATION, TAX PURPOSES, FAMILY & INDIVIDUAL - NUMBER OF OPERATIONS","Family Owned",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="GRAIN - OPERATIONS WITH SALES","Sold Grain",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="HAY & HAYLAGE - OPERATIONS WITH AREA HARVESTED","Sold Hay",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="FRUIT & TREE NUT TOTALS - OPERATIONS WITH SALES","Sold Fruit,Veges or Nuts",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="VEGETABLE TOTALS, INCL SEEDS & TRANSPLANTS, IN THE OPEN - OPERATIONS WITH SALES","Sold Fruit,Veges or Nuts",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="CATTLE, INCL CALVES - OPERATIONS WITH INVENTORY","Own Cattle",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="CHICKENS, LAYERS - OPERATIONS WITH INVENTORY","Own Chickens",short_desc)) %>%
  mutate(short_desc=ifelse(short_desc=="CHICKENS, PULLETS, REPLACEMENT - OPERATIONS WITH INVENTORY","Own Chickens",short_desc)) %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (1,000,000 TO 2,499,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (2,500,000 TO 4,999,999 $)") %>%
  filter(!domaincat_desc=="ECONOMIC CLASS: (5,000,000 OR MORE $)") %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (LESS THAN 1,000 $)","<1",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000 TO 2,499 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (2,500 TO 4,999 $)","1-5",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (5,000 TO 9,999 $)","5-10",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (10,000 TO 24,999 $)","10-25",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (25,000 TO 49,999 $)","25-50",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (50,000 TO 99,999 $)","50-100",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (100,000 TO 249,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (250,000 TO 499,999 $)","100-500",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (500,000 TO 999,999 $)","500-1,000",domaincat_desc)) %>%
  mutate(domaincat_desc=ifelse(domaincat_desc=="ECONOMIC CLASS: (1,000,000 OR MORE $)",">1,000",domaincat_desc)) %>%
  mutate(Value=as.numeric(gsub(",","",Value))/1000) %>%
  select(domaincat_desc,short_desc,Value)%>%
  group_by(domaincat_desc,short_desc) %>%
  summarise(Value=sum(Value)) %>%
  mutate(Valuelabel=ifelse(Value>50,Value,NA)) %>%
  ggplot(aes(fill=factor(domaincat_desc,levels=c("<1","1-5","5-10","10-25","25-50","50-100","100-500","500-1,000",">1,000"
  )),y=Value,x=factor(short_desc,levels=c("Sold Fruit,Veges or Nuts","Sold Hay","Sold Grain","Own Chickens","Own Cattle","Own Truck","Own Tractor","Rent Land to Others","Rent Land from Others","Own Land","Family Owned","Total")),label=floor(Valuelabel))) +
    geom_bar(position="stack",stat="identity",width=.9) +
    geom_text(size = 3, position = position_stack(vjust = .5))+
    ggtitle(paste0("Thousands of Farms by Activity"))+
    coord_flip() +
    theme_minimal()+
    scale_fill_brewer(palette="Spectral")+
    labs(x=NULL,y="Thousands of Farms",fill="Actual Sales ($000)")+
    theme(plot.title = element_text(hjust = 0.5,size=16),legend.text = element_text(hjust = 1,size=12), text = element_text(size=14))

# Draw and save plot
plot_bar_activity
ggsave("farm_activity.png")