# create a dataset
colors_in=c( "#75b9ce", "#2188a8","#0c566d")
program=c(rep("SNAP" , 3) , rep("School Lunch" , 3) , rep("WIC" , 3)  )
year=rep(c("2016" , "2017" , "2018") , 3)
value1=c(81,78,74,22,23,24,6.4,6.3,6.2)
data1=data.frame(program,year,value1)
data_test1 <- data1
data_test1$program <- as.character(data_test1$program)
data_test1$program <- factor(data_test1$program, levels = c("SNAP", "School Lunch", "WIC"))

# Grouped
ggplot(data_test1, aes(fill=year, y=value1, x=program)) + 
    geom_bar(position="dodge", stat="identity") + scale_fill_manual(values=colors_in) + 
  custom_theme + ylab("$Billions") + geom_text(aes(label=year), size=3,position=position_dodge(width=0.9), vjust=-0.25) +
   labs (title="SNAP has the most FNS funding",subtitle="FNS Budget Authority (billions of dollars) 2016 - 2018",  caption="Source:USDA Food Nutrition Services") + theme(legend.position="none")

value2=c(44.2,42.8,42.2,30.4,30.6,30.6,7.7,7.3,7.2)
data2=data.frame(program,condition,value2)

data_test2 <- data2
data_test2$program <- as.character(data_test2$program)
data_test2$program <- factor(data_test2$program, levels = c("SNAP", "School Lunch", "WIC"))

ggplot(data_test2, aes(fill=year, y=value2, x=program)) + ylab("People Served(millions)") + 
    geom_bar(position="dodge", stat="identity") + scale_fill_manual(values=colors_in) + custom_theme + labs (title="SNAP served more people than school lunch and WIC combined",subtitle="People Served Through Nutrition Assistance Programs 2016 - 2018",caption="Source:USDA Food Nutrition Services") + geom_text(aes(label=year), size=3,position=position_dodge(width=0.9), vjust=-0.25) + theme(legend.position="none")

The above two graphs show the budget distribution for FNS (Food nutrition services) and number of people served by those three nutrition programs. In specific, SNAP takes three times more FNS budget than WIC and School Lunch program. And the number of participants of SNAP are roughly the WIC and School Lunch combined. Across the three years from 2016 to 2018, the budget to the WIC and SNAP declined, and correspondingly, the number of participants declined as well. A little bit more funds went to School Lunch program and the participants increased over the year.

snap_event<-read.xlsx("snap_event.xlsx",head=T,sheetIndex = 1)

snap_event$y<-c(runif(8,.3,.5),runif(7,-.5,-.3))

positions <- c(.5, -.5, 1.0, -1.0, 1.5, -1.5)
directions <- c(1, -1)

line_pos <- data.frame(
    "year"=unique(snap_event$year),
    "position"=rep(positions, length.out=length(unique(snap_event$year))),
    "direction"=rep(directions, length.out=length(unique(snap_event$year)))
)

snap_event <- merge(x=snap_event, y=line_pos, by="year", all = TRUE)
snap_event <- snap_event[with(snap_event, order(year, event)), ]

text_offset <- 0.05
snap_event$text_position <- (text_offset * snap_event$direction) + snap_event$position

ggplot(snap_event,aes(x=year,y=0,  label=event))+labs(col="events")+theme_classic()+geom_hline(yintercept=0, 
                color = "black", size=1.5)+geom_segment(data=snap_event, aes(y=position,yend=0,xend=year), color='black', size=0.6)+
scale_y_discrete(limits=c(-1,1))+theme(axis.line.y=element_blank(),
                 axis.text.y=element_blank(),
                 axis.title.x=element_blank(),
                 axis.title.y=element_blank(),
                 axis.ticks.y=element_blank(),
                 #axis.text.x =element_blank(),
                 axis.ticks.x =element_blank(),
                 axis.line.x =element_blank()) + scale_x_continuous(breaks =c(1970,1980,1990,2000,2010,2020)) + 
  # geom_text(data=snap_event, aes(x=year,y=-0.2,label=year),vjust=1.4,size=1.5, color="#2188a8")+
  
  geom_point(aes(y=0), size=2,color="#2188a8")+
  
  geom_text(aes(y=text_position,label=event),size=1.9,fontface="bold",vjust=0)+labs(title="Big Events in SNAP History", subtitle="SNAP 1964 - Today", caption="Source:USDA Food Nutrition Services") + theme(
          plot.title = element_text(size = 11,face = "bold",
                                  family = "Garamond",
                              lineheight = 0.5,
                                   color = textcolor),
       plot.subtitle = element_text(size = 11,
                                  family = "Garamond",face = "italic",
                              lineheight = 0.5,
                                   color = textcolor),
        plot.caption = element_text(size = 6,
                                  family = "Arial",
                              lineheight = 0.5,
                                   color = textcolor),
         plot.margin =          margin(t = 30, 
                                       r = 30, 
                                       b = 30, 
                                       l = 30, 
                                    unit = "pt"),
          axis.title = element_text(size = 13,
                                    face = "bold",
                                  family = "Arial",
                              lineheight = 0.5,
                                   color = textcolor),
           axis.text = element_text(size = 5,
                                    face = "bold",
                                  family = "Arial",
                              lineheight = 0.5, 
                                   color = textcolor))

From 1989 to 2018, the SNAP participation number in US increased from below two millions to about four millions, almost doubled in two decades. The peaks within the two decades happened at year 1994 and year 2013, reached over 2.5 millions and over 4.5 millions respectively. Connecting to the previous SNAP timeline graph, SNAP reached participation milestone for those two years. The Food Stamp Nutrition Education program at 1992 and Mickey Leland Childhood Hunger Relief Act at 1993 influenced the growth of SNAP, the participants increased dramatically afterwards. And from 2002 to 2014, there are five events happened that drove the SNAP participaion greatly.

# skip first 6 rows, and only read in 1st and 4th col
part_18 <-read.xlsx("snap_part_2018.xls",header=F,skip=6,sheetIndex = 1)
part_18<-part_18[-c(9,48,49),][7:58,]
colnames(part_18) <- c("state","sep17","aug18","sep18","sep18_vs_aug18","sep18_vs_sep17")
r1="Northeast"
r2="Midwest"
r3="South"
r4="West"

s1<-c("Connecticut", "Maine", "Massachusetts", "New Hampshire", "Rhode Island", "Vermont", "New Jersey", "New York","Pennsylvania")
s2<-c("Illinois", "Indiana", "Michigan", "Ohio", "Wisconsin", "Iowa", "Kansas", "Minnesota", "Missouri", "Nebraska", "North Dakota", "South Dakota")
s3<-c("Delaware", "Florida", "Georgia", "Maryland", "North Carolina", "South Carolina", "Virginia", "District of Columbia", "Alabama", "Kentucky", "Mississippi", "Tennessee", "Arkansas", "Louisiana", "Oklahoma", "Texas")
s4<-c("Arizona", "Colorado", "Idaho", "Montana", "Nevada", "New Mexico", "Utah", "Wyoming","Alaska", "California", "Hawaii", "Oregon", "Washington")


for (row in 1:nrow(part_18)) {
  state <- part_18[row, "state"]
  if (state %in% s1) {
  region = r1
}else if(state %in% s2){
  region = r2
}else if(state %in% s3){
  region = r3
}else if(state %in% s4) {
  region = r4
}
  part_18$region[row] <- region}
data = data %>% arrange(region, as.numeric(as.character(sep18)))

# Set a number of 'empty bar' to add at the end of each group
empty_bar=3
to_add = data.frame( matrix(NA, empty_bar*nlevels(data$region), ncol(data)) )
colnames(to_add) = colnames(data)
to_add$group=rep(levels(data$region), each=empty_bar)
data=rbind(data, to_add)
data=data %>% arrange(region)
data$id=seq(1, nrow(data))


# Get the name and the y position of each label
label_data=data
number_of_bar=nrow(label_data)
angle= 90 - 360 * (label_data$id-0.5) /number_of_bar     # I substract 0.5 because the letter must have the angle of the center of the bars. Not extreme right(1) or extreme left (0)
label_data$hjust<-ifelse( angle < -90, 1, 0)
label_data$angle<-ifelse(angle < -90, angle+180, angle)
 
# prepare a data frame for base lines
base_data=data %>% 
  group_by(region) %>% 
  summarize(start=min(id), end=max(id) - empty_bar) %>% 
  rowwise() %>% 
  mutate(title=mean(c(start, end)))


# prepare a data frame for grid (scales)
grid_data = base_data
grid_data$end = grid_data$end[ c( nrow(grid_data), 1:nrow(grid_data)-1)] + 1
grid_data$start = grid_data$start - 1
grid_data=grid_data[-1,]

colors_in = c("#2188a8","#dd8944","#11af80","#db6bac")

# Make the plot
p = ggplot(data, aes(x=as.factor(id), y=as.numeric(as.character(sep18)), fill=region))+ scale_fill_manual(values=colors_in)+       # Note that id is a factor. If x is numeric, there is some space between the first bar

  geom_bar(aes(x=as.factor(id), y=as.numeric(as.character(sep18)), fill=region), stat="identity", alpha=0.5) +

  # Add a val=100/75/50/25 lines. I do it at the beginning to make sur barplots are OVER it.
  geom_segment(data=grid_data, aes(x = end, y = 4000000, xend = start, yend = 4000000), colour = "grey", alpha=1, size=0.1 , inherit.aes = FALSE ) +
  geom_segment(data=grid_data, aes(x = end, y = 3000000, xend = start, yend = 3000000), colour = "grey", alpha=1, size=0.1 , inherit.aes = FALSE ) +
  geom_segment(data=grid_data, aes(x = end, y = 2000000, xend = start, yend = 2000000), colour = "grey", alpha=1, size=0.1 , inherit.aes = FALSE ) +
  geom_segment(data=grid_data, aes(x = end, y = 1000000, xend = start, yend = 1000000), colour = "grey", alpha=1, size=0.1 , inherit.aes = FALSE ) +

  # Add text showing the value of each 100/75/50/25 lines
  annotate("text", x = rep(max(data$id),4), y = c(1000000, 2000000, 3000000, 4000000), label = c("1m", "2m", "3m", "4m") , color="grey", size=3 , angle=0, fontface="bold", hjust=5) +
  geom_bar(aes(x=as.factor(id), y=as.numeric(as.character(sep18)), fill=region), stat="identity", alpha=0.5) +
  ylim(-1500000,4000000) +
theme_minimal() +
theme(plot.title = element_text(size = 11,face = "bold",
                                  family = "Garamond",
                              lineheight = 0.5,
                                   color = textcolor),
       plot.subtitle = element_text(size = 11,
                                  family = "Garamond",face = "italic",
                              lineheight = 0.5,
                                   color = textcolor),
        plot.caption = element_text(size = 6,
                                  family = "Arial",
                              lineheight = 0.5,
                                   color = textcolor),
                              legend.position = "none",
      axis.text = element_blank(),
   axis.title = element_blank(),
   panel.grid = element_blank(), plot.margin = margin(t = 30,
                                       r = 30,
                                       b = 30,
                                       l = 10,
                                    unit = "pt"))+

  coord_polar()+
  geom_text(data=label_data, aes(x=id, y=as.numeric(as.character(sep18))+100000, label=state, hjust=hjust), color="black", fontface="bold",alpha=0.6, size=2, angle= label_data$angle, inherit.aes = FALSE ) +
  # Add base line information
  geom_segment(data=base_data, aes(x = start, y =-100000, xend = end, yend = -100000), colour = "black", alpha=0.8, size=0.6, inherit.aes = FALSE )  +
  geom_text(data=base_data, aes(x = title, y = -180000, label=region), hjust=c(0.9,1,0.35,0), colour = "black", alpha=0.8, size=2, fontface="bold", inherit.aes = FALSE)+
  labs(title="SNAP participation ranges widely across states in different regions", subtitle="US SNAP Participation population by State in 2018" ) 

p3<-add_sub(p, "Source: USDA Food Nutrition", y  = 0, hjust=-1.2,vjust = -5,size=6.5) 
ggdraw(p3)

# annotate("text",x=-Inf,y=-Inf,hjust=5,vjust=0,label="Source: USDA Food Nutrition")
# caption = "Source: USDA Food Nutrition Services\nhttps://www.fns.usda.gov/pd/supplemental-nutrition-assistance-program-snap"

Participation population of SNAP program in US across the 51 states, four regions, are shown above. West has the widest range of participation numbers, while the it’s more concentrated for Midwest. Overall, at state level, California has the most participations, nearly 4 million, and Wyoming has the fewest participants, about 0.1 million.

library(reshape2)
meltdf <- melt(time.series,id="Year")

ggplot(meltdf,aes(x=Year,y=as.numeric(value),colour=variable,group=variable)) + ylab("% Change") + geom_line() +geom_point() +
  labs(title = "SNAP participation change shares similar trend with poverty rate",
        subtitle = "SNAP Participants and Economic Indicators: Percentage Change 1969 - 2018",
         caption = "Source:Bureau of Labor Statistics"
    ) +
  scale_color_manual(values = c("#dd8944","#2188a8","#db6bac"))+scale_x_discrete(breaks =c(1970,1980,1990,2000,2010,2020)) +
scale_y_continuous(limits=c(-.2,1.2)) + 
  custom_theme 

Overall, the changes of SNAP participation population, unemployment rate, and poverty rate have similar trends. Unemployment rate change fluctuates more than the other two rates, while the poverty rate change is the most consistent over the years. One outlier is in 1971, the SNAP population doubled compared to the previous year and then came back.

um<-read.xlsx("um_79to18.xlsx",header=TRUE,sheetIndex = 1)[11:40,]
part_ben$um<-um$Unemployment.Rate

b <- ggplot(part_ben, aes(x = Year, y=Benefit.dollar.person.))

b + geom_point(aes(colour = um), size = 5)+ scale_color_gradientn(colours = c("#a5bfc6","#2188a8","#0c566d","#033d4c","#dd8944"))+ scale_x_continuous(breaks =c(1990,1995,2000,2005,2010,2015))+
  
  #scale_fill_gradientn(colours = c("#a5bfc6","#2188a8","#0c566d","#033d4c","#dd8944"), breaks=c(4,6,8,10))+ 
  
  geom_vline(xintercept = 2008, color = "red", linetype = "dashed") +
  annotate("text", label = "Recession", x = 2008, y = 110, family="Garamond",color = "black")+ labs(y = "Benefits per person ($)",x = "Year",subtitle="SNAP Benefits and Unemployment Rate 1989-2018", title="SNAP participation and benefits increased over the last two decades",caption="Source:Bureau of Labor Statistics") +custom_theme+theme(
            legend.text = element_text(size = 8),
 legend.title=element_text(size=6,face="bold")) 

#par(mar=c(3, 3, 3, 3))

From 1989 to 2018, the benefits per person increased dramatically. From around 50 dollars to above 250 dollars, the dollar benefitfive is currenly five times more than from the beginning. The color scale shows the range of unemployment rate over the years. The unemployment rate fluctuates from 1989 to 2008, and increased dramatically around 2008 to 2014, which is the financial crisis period. The benefit per person decreased a little, but after the crisis, the benefit jumped very high till around 250 dollars per person and came back to around 125 in the most recent year.

pop <- read.delim("data.txt",header = TRUE, sep = ",")
# part_18:
# Guam
# Virgin Islands

# pop:
# Rhode Island

# pop$part.18
# pop$X2019.Population
# pop$part.per.capita

colnames(pop)[colnames(pop)=="State"] <- "state_name"
# part_18<-part_18[ ! part_18$state %in% c("Guam","Virgin Islands"), ]
# pop<              -pop[ ! pop$state %in% c("Rhode Island"), ]
# pop$part.18<-part_18$sep18
# pop$part.per.capita<-as.numeric(pop$part.18)/as.numeric(pop$X2019.Population)

pop %>%
  left_join(states, by = "state_name") %>% 
  ggplot(mapping = aes(long, lat, group = group, fill = X..of.US)) +labs(subtitle="US SNAP Participation per capita by State 2018", title="California, Texas and Florida have the highest SNAP participation rate", caption="Source:USDA Food Nutrition Services")+custom_theme+theme(legend.text = element_text(size = 8),legend.title=element_text(size=7,face="bold"),panel.grid.major = element_line(colour = 'transparent'),
 axis.title.x=element_blank(), 
 axis.text.x=element_blank(),
 axis.ticks.x=element_blank(),
 axis.title.y=element_blank(), 
 axis.text.y=element_blank(),
 axis.ticks.y=element_blank(),
 panel.background=element_blank(),
 panel.border=element_blank(),
 panel.grid.minor=element_blank(),
 plot.background=element_blank(),
 legend.position = c(.95, .5)
 )+
  geom_polygon(color = "black", size = 0.25) + scale_fill_gradientn(colours = c("#a5bfc6","#2188a8","#0c566d","#033d4c","#dd8944"), breaks=c(0,.02,.04,.06,.08,.1,.12)) + 
                                                                    
  coord_map("albers", lat0 = 39, lat1 = 45)  + labs(fill = "Participation per capita")

The first map is the US base map with state names labeled. From the US state map, combining with SNAP participation data, we can see that state AL, ND, SD and NH have the highest participation rate. This is reasonable as these states are less developed and tend to have more poor people. The participation rate range is not very big, therefore the colors are generally darker. Most of the states have participation rate around 0.00002 or lower.

dt %>%
  left_join(states, by = "state_abbv") %>% 
  ggplot(mapping = aes(long, lat, group = group, fill = FOOD_TAX14))  +labs(subtitle="US Food Tax 2018", title="Most of the states have low food tax", caption="Source:USDA Food Nutrition Services
")+custom_theme+theme(legend.text = element_text(size = 8),legend.title=element_text(size=7,face="bold"),panel.grid.major = element_line(colour = 'transparent'),
      
 axis.title.x=element_blank(), 
 axis.text.x=element_blank(),
 axis.ticks.x=element_blank(),
 axis.title.y=element_blank(), 
 axis.text.y=element_blank(),
 axis.ticks.y=element_blank(),
 panel.background=element_blank(),
 panel.border=element_blank(),
 panel.grid.minor=element_blank(),
 plot.background=element_blank(),
 legend.position = c(.9, .5)
 )+
  geom_polygon(color = "black", size = 0.25) + scale_fill_gradientn(colours = c("#a5bfc6","#2188a8","#0c566d","#033d4c","#dd8944"), breaks=c(0,2,4,6,8)) + 
  coord_map("albers", lat0 = 39, lat1 = 45)  + labs(fill = "Food Tax") 

Drawing a US Food Tax map across states and comparing it with the previous participation rate map, we can see that places with higher SNAP participation rates do not have high food tax. And places with high food tax rates tend to have low SNAP participation rates.

# head(retailer)
cnt <- retailer %>% group_by(State) %>% count(State)

colnames(cnt)<-c("state","store_count")
statebins_continuous(cnt, value_col="store_count",brewer_pal="Blues", font_size = 3,
                     legend_title="Store counts")+
labs(title="SNAP Retail Store Counts by State", caption="Source:USDA Food Nutrition Services
")+custom_theme+
  theme(legend.text = element_text(size = 6),legend.title=element_text(size=7,face="bold"),panel.grid.major = element_line(colour = 'transparent'),
 axis.title.x=element_blank(), 
 axis.text.x=element_blank(),
 axis.ticks.x=element_blank(),
 axis.title.y=element_blank(), 
 axis.text.y=element_blank(),
 axis.ticks.y=element_blank(),
 panel.background=element_blank(),
 panel.border=element_blank(),
 panel.grid.minor=element_blank(),
 plot.background=element_blank()) + scale_fill_gradientn(colours = c("#a5bfc6","#dd8944"), breaks=c(0,5000,10000,15000,20000))

This above graph aggregates the county retailer numbers to state level. And as shown above, we can see that CA has the most number of retailer stores, followed by NY, TX and FL. The most number of stores in one state reaches about 25000, while the least number is about 1500. The states with few number of stores mostly locate at north west of the states.

snap_policy_index<-read.xlsx("snappolicyindexdata.xls",sheetIndex = 1,header=TRUE)[2:971,][c(2,4,10)]
colnames(snap_policy_index)<-c("state_abbv","year","weighted_policy_index")
# head(snap_policy_index)
dt96<-snap_policy_index %>% group_by("state_abbv") %>% filter(year=="1996")
dt14<-snap_policy_index %>% group_by("state_abbv") %>% filter(year=="2014")
dt<- merge(dt96, dt14, by="state_abbv", all=TRUE)

dt$index_diff<-as.numeric(as.character(dt$weighted_policy_index.y))-as.numeric(as.character(dt$weighted_policy_index.x))
# dt
# states
dt %>%
  left_join(states, by = "state_abbv") %>%
ggplot(mapping = aes(long, lat, group = group, fill = index_diff))  +labs(title="SNAP policies became more accommodative toward enrollment in all\nthe states, North Carolina and New York the most", subtitle="SNAP policy index changes 1996 - 2014", caption="Source:USDA Food Nutrition Services
")+custom_theme+theme(panel.grid.major = element_line(colour = 'transparent'),
 axis.title.x=element_blank(), 
 axis.text.x=element_blank(),
 axis.ticks.x=element_blank(),
 axis.title.y=element_blank(), 
 axis.text.y=element_blank(),
 axis.ticks.y=element_blank(),
 panel.background=element_blank(),
 panel.border=element_blank(),
 panel.grid.minor=element_blank(),
 plot.background=element_blank(),
 legend.position = c(.9, .5),
 legend.text = element_text(size = 8),
 legend.title=element_text(size=7,face="bold")
 )+
  geom_polygon(color = "black", size = 0.25) + scale_fill_gradientn(colours = c("#a5bfc6","#2188a8","#0c566d","#033d4c","#dd8944"), breaks=c(2,2.5,3,3.5,4)) + 
  coord_map("albers", lat0 = 39, lat1 = 45)  + labs(fill = "Index Differences") 

The graph above depicts changes in the SNAP Policy Index across the States from 1996 to 2014. Index values range from 1 to 10, with 10 representing a set of policies that is the most accommodative to SNAP participation. The smallest change occurred in Wyoming, where the index went from a value of around 4.1 to 6.6. The largest change occurred in New York, where the index moved from near 2.3 to over 8.6. It is notable that in all 50 States and the District of Columbia, the index increased over the period from 1996 to 2014, meaning that SNAP policies became more accommodative toward enrollment in SNAP in every one of those States.

# part_ben$Year <- as.Date(as.character(part_ben$Year), format = "%Y")
# snap_event$year <- as.Date(as.character(snap_event$year), format = "%Y")
# 
# yrng <- range(part_ben$Participation.number.)
# ymin <- yrng[1]
# ymax <- yrng[1] + 0.1*(yrng[2]-yrng[1])
# 
# p2 <- ggplot()
# p2 <- p2 + geom_line(mapping=aes(x=Year, y=Participation.number.), data=part_ben , size=2, alpha=0.5,color='blue') 
# p2 <- p2 + scale_x_date("time")  +
#   # scale_x_discrete(limits=c("1990","1995","2000","2005","2010","2015"))+
# scale_y_continuous(name="participation")+ labs(title="SNAP Participation Population 1989-2018", subtitle="SNAP participation doubled over the last two decades", caption="Source:USDA Food Nutrition Services\nhttps://www.fns.usda.gov/snap/short-history-snap")+theme(plot.caption = element_text(hjust = 0),plot.title=element_text(size=20, 
#                                     family="Garamond",
#                                     color="black",
#                                     hjust=0.5,
#                                     lineheight=1.2),  # title
#             plot.subtitle=element_text(size=15, 
#                                        family="Garamond",
#                                        face="italic",
#                                        hjust=0.5))
# p2

FNS Budget provides funding for the major nutrition assistance programs, accounting for projected program participation and food cost inflation. It seeks to prevent and reduce food insecurity and to improve the nutritional status of recipients.

To support FNS’ work to identify and eliminate fraud, waste, and abuse, the Budget provides resources for program integrity efforts in all of the major programs, including SNAP, WIC, and the Child Nutrition Programs.

data=as.data.frame(matrix(c(10.3,11.7,19.9,24,32.1,31,26.4,24.8,11.3,8.5), ncol=5))
# data
colnames(data)=c("Excellent", "Very Good" , "Good" , "Fair" , "Poor" )
rownames(data)=c("Begining","AfterSixMonth")

# To use the fmsb package, I have to add 2 lines to the dataframe: the max and min of each topic to show on the plot!
data=rbind(rep(35,5) , rep(0,5) , data)
colors_border=c( rgb(0.2,0.5,0.5,0.9), rgb(0.8,0.2,0.5,0.9) , rgb(0.7,0.5,0.1,0.9))
colors_in=c( rgb(0.2,0.5,0.5,0.4), rgb(0.8,0.2,0.5,0.4) , rgb(0.7,0.5,0.1,0.4) )

radarchart(data , axistype=1,
    #custom polygon
    pcol=colors_border , pfcol=colors_in , plwd=2 , plty=1,family = "Arial",
    #custom the grid
    cglcol="grey", cglty=1, axislabcol="grey", caxislabels=seq(0,40,8), cglwd=1,
    #custom labels
    vlcex=0.8 
    )
legend(x=0.7, y=1, legend = rownames(data[-c(1,2),]), bty = "n", pch=20 , col = colors_in, text.col = "black", cex=1, pt.cex=1)
title(main=expression(paste(bold("Health Status improved after 6 months participation in SNAP"))),family="Garamond", font.main= 1,cex.main=1.2)
text(-4,-1, "Source:Bureau of Labor Statistics",cex = .8,
     adj = c(-2.5,3))

par(mar=c(3, 3, 3, 3))
#title(main=expression(paste(bold("Hello"), " world")), font.main=1)

# add_sub(p, "Source:Bureau of Labor Statistics", y  = 0, hjust=-.4,vjust = -12) 
# p3<-add_sub(p, "Source:Bureau of Labor Statistics", y  = 0, hjust=-.4,vjust = -12) 
# ggdraw(p3)
#fig.cap="Source:Bureau of Labor Statistics"

The longitudinal estimates compare new SNAP participants to the same participants about six months later. Tabulations are based on the following overall sample sizes: 3,275 new-entrant households and 3,375 six-month households in the cross-sectional sample; and 3,275 new-entrant households observed at baseline and again at follow up six months later in the longitudinal sample. We can see that the orange area leans towards good side more than the blue area, which suggests that the health status after six month for SNAP participants is better than the start of the survey.