#"ACEs_scoping_review"

# PRELIMINARIES 

## This chunk loads the data file, creates various sub-data frames & variables, & runs a few home-made functions 

#required packages
library(ggplot2)
library(tableone)
library(gmodels)
library(rworldmap)

#First set the working directory to the location of the "Lit_reviewR.csv" file
setwd("C:/Users/karat004/Documents/Work stuff/Child Welfare Stuff/Lit review")

#to import & create data sets & variables
ScopingR<- read.csv("Lit_reviewR.csv") #copy and paste as values all the data from the "To be exported" sheet of the Excel file to a csv file and use your own path and file name here
Scoping_include<-subset(ScopingR,Exclude_or_include=="INCLUDE") #includes only the INCLUDED articles
Scoping_primary<-subset(Scoping_include, Purpose_secondary!="yes") #excludes the EXCLUDED articles & all secondary-purpose articles

Scoping_1999<-subset(Scoping_include,Year=="1999")
Scoping_2000<-subset(Scoping_include,Year=="2000")
Scoping_2001<-subset(Scoping_include,Year=="2001")
Scoping_2002<-subset(Scoping_include,Year=="2002")
Scoping_2003<-subset(Scoping_include,Year=="2003")
Scoping_2004<-subset(Scoping_include,Year=="2004")
Scoping_2005<-subset(Scoping_include,Year=="2005")
Scoping_2006<-subset(Scoping_include,Year=="2006")
Scoping_2007<-subset(Scoping_include,Year=="2007")
Scoping_2008<-subset(Scoping_include,Year=="2008")
Scoping_2009<-subset(Scoping_include,Year=="2009")
Scoping_2010<-subset(Scoping_include,Year=="2010")
Scoping_2011<-subset(Scoping_include,Year=="2011")
Scoping_2012<-subset(Scoping_include,Year=="2012")
Scoping_2013<-subset(Scoping_include,Year=="2013")
Scoping_2014<-subset(Scoping_include,Year=="2014")
Scoping_2015<-subset(Scoping_include,Year=="2015")
Scoping_2016<-subset(Scoping_include,Year=="2016")
Scoping_2017<-subset(Scoping_include,Year=="2017")
Scoping_2018<-subset(Scoping_include,Year=="2018")
Scoping_2019<-subset(Scoping_include,Year=="2019")
Scoping_2020<-subset(Scoping_include,Year=="2020")

Scoping_primary_1999<-subset(Scoping_primary,Year=="1999")
Scoping_primary_2000<-subset(Scoping_primary,Year=="2000")
Scoping_primary_2001<-subset(Scoping_primary,Year=="2001")
Scoping_primary_2002<-subset(Scoping_primary,Year=="2002")
Scoping_primary_2003<-subset(Scoping_primary,Year=="2003")
Scoping_primary_2004<-subset(Scoping_primary,Year=="2004")
Scoping_primary_2005<-subset(Scoping_primary,Year=="2005")
Scoping_primary_2006<-subset(Scoping_primary,Year=="2006")
Scoping_primary_2007<-subset(Scoping_primary,Year=="2007")
Scoping_primary_2008<-subset(Scoping_primary,Year=="2008")
Scoping_primary_2009<-subset(Scoping_primary,Year=="2009")
Scoping_primary_2010<-subset(Scoping_primary,Year=="2010")
Scoping_primary_2011<-subset(Scoping_primary,Year=="2011")
Scoping_primary_2012<-subset(Scoping_primary,Year=="2012")
Scoping_primary_2013<-subset(Scoping_primary,Year=="2013")
Scoping_primary_2014<-subset(Scoping_primary,Year=="2014")
Scoping_primary_2015<-subset(Scoping_primary,Year=="2015")
Scoping_primary_2016<-subset(Scoping_primary,Year=="2016")
Scoping_primary_2017<-subset(Scoping_primary,Year=="2017")
Scoping_primary_2018<-subset(Scoping_primary,Year=="2018")
Scoping_primary_2019<-subset(Scoping_primary,Year=="2019")
Scoping_primary_2020<-subset(Scoping_primary,Year=="2020")

# this is a function for calculating proportion of outliers
outlierSummary<-function(variable, digits = 2){
  zvariable<-(variable-mean(variable, na.rm = TRUE))/sd(variable, na.rm = TRUE)
  outlier95<-abs(zvariable) >= 1.96
  outlier99<-abs(zvariable) >= 2.58
  outlier999<-abs(zvariable) >= 3.29
  ncases<-length(na.omit(zvariable))
  percent95<-round(100*length(subset(outlier95, outlier95 == TRUE))/ncases, digits)
  percent99<-round(100*length(subset(outlier99, outlier99 == TRUE))/ncases, digits)
  percent999<-round(100*length(subset(outlier999, outlier999 == TRUE))/ncases, digits)
  cat("Absolute z-score greater than 1.96 = ", percent95, "%", "\n")
  cat("Absolute z-score greater than 2.58 = ", percent99, "%", "\n")
  cat("Absolute z-score greater than 3.29 = ", percent999, "%", "\n")
}

#this is a function for generating frequency & proportion data

mytable<-function(x, exclude="") {
  freq<-table(x, exclude="")
  prop<-round(prop.table(table(x, exclude=""))*100,1)
  cbind(Freq=freq, Percent=prop)
}


#this is a function for generating my favorite stats
mystats<-function(x,...){
  N<- length(!is.na(x))
  Mean<-ifelse(sum(!is.na(x))==0, NA, round(mean(x, na.rm=TRUE),2))
  SD<-ifelse(sum(!is.na(x))<=1, NA, round(sd(x,na.rm=TRUE),2))
  Median<-ifelse(sum(!is.na(x))==0, NA, round(median(x, na.rm=TRUE),2))
  Min<-ifelse(sum(!is.na(x))==0, NA, round(min(x,na.rm=TRUE),2))
  Max<-ifelse(sum(!is.na(x))==0, NA, round(max(x,na.rm=TRUE),2))
  c(N = N, Mean = Mean, SD = SD, Median = Median, Min = Min, Max = Max)
}


#this is a function for just generating mean & sd
mean.and.sd <- function(x=1:10){ 
  av <- round(mean(x, na.rm=TRUE),1)
  sd <- round(sqrt(var(x,na.rm=TRUE)),1)
  c(mean=av, SD=sd) 
}


#this is a function for doing exploratory data analysis w/continous variables
explore_cont<-function(x){
  oldpar<-par(no.readonly=TRUE,mfrow=c(2,2),mar=c(2,4,1,2))
  hist(x, main="")
  plot(density(x, na.rm=TRUE), xlim=c(min(x, na.rm=TRUE),max(x,na.rm=TRUE)), type = "l", main="")
  qqnorm(y=x, main="");qqline(y=x, distribution=qnorm)
  boxplot(x,main="")
  print(str(x))
  print(round(stat.desc(x,basic=TRUE, norm=TRUE),1))
  outlierSummary(x)
  par(oldpar)
}

#this is a function for doing exploratory data analysis w/categorical variables
explore_categ<-function(x) {
  print(attributes(x))
  print(mytable(x))
  cat("\n")
  cat("Total N = ", length(x))
}


# This is a function to add numbers to boxplots
give.n <- function(x){
  return(c(y = median (x)*2.5, label = length(x))) 
  # goes in this command: stat_summary(fun.data = give.n, geom = "text")+# min 0.7
}

## This chunk checks the data file to make sure everything looks in order 

class(ScopingR)
dim(ScopingR)
str(ScopingR)
names(ScopingR)

# of included & excluded & primary/secondary articles
explore_categ(ScopingR$Exclude_or_include) # of included & excluded articles
explore_categ(Scoping_include$Purpose_secondary) # of secondary articles
1361-80 # Number of primary purpose articles = 1281


#Year of publication

explore_categ(Scoping_include$Year) 
Scoping_include$Year_f<-factor(Scoping_include$Year, levels = c(1999,2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019, 2020), labels = c("1999", "2000", "2001", "2002", "2003", "2004", "2005", "2006", "2007", "2008", "2009", "2010", "2011", "2012", "2013", "2014", "2015", "2016", "2017","2018", "2019", "2020")) 
explore_categ(Scoping_include$Year_f) 

377/21 #18-fold increase from 2009 to 2019

#Journals in which the articles were published

length(unique(Scoping_include$Journal)) # of unique journals
explore_categ(Scoping_include$Journal) # number and proportion of articles in each journal
13 + 16 + 22 + 16 + 17 + 23 + 17 + 26 + 15 + 18 + 17 + 19 + 13 # number of articles in the 17 journals (other than Child Abuse & Neglect) with at least 1% of the articles
1 + 1.2 + 1.6 + 1.2 + 1.3 + 1.7 + 1.3 + 1.9 + 1.1 + 1.3 + 1.3 + 1.4 + 1.0 # proportion of articles in the 17 journals (other than Child Abuse & Neglect) with at least 1% of the articles


length(unique(Scoping_1999$Journal))
length(unique(Scoping_2000$Journal))
length(unique(Scoping_2001$Journal))
length(unique(Scoping_2002$Journal))
length(unique(Scoping_2003$Journal))
length(unique(Scoping_2004$Journal))
length(unique(Scoping_2005$Journal))
length(unique(Scoping_2006$Journal))
length(unique(Scoping_2007$Journal))
length(unique(Scoping_2008$Journal))
length(unique(Scoping_2009$Journal))
length(unique(Scoping_2010$Journal))
length(unique(Scoping_2011$Journal))
length(unique(Scoping_2012$Journal))
length(unique(Scoping_2013$Journal))
length(unique(Scoping_2014$Journal))
length(unique(Scoping_2015$Journal))
length(unique(Scoping_2016$Journal))
length(unique(Scoping_2017$Journal))
length(unique(Scoping_2018$Journal))
length(unique(Scoping_2019$Journal))
length(unique(Scoping_2020$Journal))


(Figure2a<-ggplot(Scoping_include, aes(x = Unique_year))+
    geom_line (aes(y = No_articles), size = 1.5, linetype = "longdash", color = "black")+
    geom_line (aes(y = Unique_journal), size = 1.5, linetype = "dotted", color = "black", position = "jitter")+
    scale_y_continuous(name = "Number of articles/year", limits= c(0, 400), breaks = seq(from = 0, to = 400, by = 20),sec.axis = sec_axis(~.*1,breaks = seq(from = 0, to = 400, by = 20), name = "Number of journals/year"))+
    scale_x_continuous(limits= c(1999, 2019), breaks = seq(from = 1999, to = 2019, by = 2))+
    theme(axis.title.y = element_text(color = "black", size=18, vjust=5),axis.title.y.right = element_text(color = "black", size=16,vjust=5)) +
    theme(axis.text.y = element_text(color="black",size=14))+
    theme(axis.title.x = element_text (color = "black",  size = 18,vjust=-4),axis.text.x=element_text(size=14))+
    labs(x = "Date of publication")+
    theme(panel.grid.major = element_line(colour="gray",linetype="dotted", size=0.1),panel.background = element_rect(fill="white") ,panel.border = element_rect(color = "black", fill = NA, size =1))+
    theme(plot.background=element_rect(fill="white"),plot.margin = unit(c(2, 1, 2, 1), "cm"))) 

ggsave("Figure2a.jpeg", dpi=800, height=8, width=9.5, units="in")       

#Study Purpose

Purpose_for_table<-c("Purpose_descriptive_relational", "Purpose_screening","Purpose_intervention", "Purpose_methodological","Purpose_secondary")
Purpose_table<-CreateCatTable(vars=Purpose_for_table,data=Scoping_include, test=F)
print(Purpose_table, showAllLevels = FALSE)

Purpose_table_by_year_f<-CreateCatTable(vars=Purpose_for_table,strata = c("Year_f"), data=Scoping_include, test=F)
print(Purpose_table_by_year_f, strata = c("Year_f"), showAllLevels = FALSE)

Purpose_data<- subset(ScopingR, ScopingR$Purpose_year != "NA" & ScopingR$Purpose_category != "NA" & ScopingR$Purpose_number != "NA" & ScopingR$Purpose_proportion_year != "NA", select=c("Purpose_year", "Purpose_category", "Purpose_number", "Purpose_proportion_year"))

(Figure3<-ggplot(Purpose_data, aes(x = Purpose_year, y = Purpose_proportion_year, fill = Purpose_category, label = Purpose_number)) +
    geom_bar(stat = "identity") +
    geom_text(size = 5.5, position = position_stack(vjust = 0.5))+
    labs(x = "Year of publication", y = "Proportion of articles/year")+
    scale_x_continuous(breaks=seq(1999,2019,2))+
    scale_y_continuous(breaks=seq(0, 1.0, by=0.1))+
    scale_fill_manual(values = c("Descriptive/relational" = "gray95", "Intervention" = "gray75", "Methodological" = "gray60", "Screening"= "gray45", "Secondary" = "gray35"))+
    guides(fill=guide_legend(title="Study purpose"))+
    theme(legend.title = element_text (size=18))+
    theme(legend.text=element_text(size=16))+
    theme(axis.text=element_text(color="black",size=16))+
    theme(axis.title.x = element_text (color = "black",  size = 20,vjust=-4))+
    theme(axis.title.y = element_text(color = "black", size=20, vjust=5)) +
    theme(panel.grid.major = element_line(colour="white",linetype="dashed", size=0.1),panel.background = element_rect(fill="white") ,panel.border = element_rect(color = "black", fill = NA)))+
  theme(plot.background=element_rect(fill="white"),plot.margin = unit(c(2, 1, 2, 1), "cm")) 

ggsave("Figure3.jpeg", dpi=800, height=10, width=14, units="in")      


#Operational definitions of ACEs

ACEs_for_table<-c("ACE_emotional_abuse", "ACE_emotional_neglect", "ACE_physical_abuse","ACE_sexual_abuse", "ACE_physical_neglect", "ACE_other_maltreatment",
                  "ACE_mental_illness","ACE_substance","ACE_legal","ACE_domestic_violence","ACE_divorce","ACE_parental_instability","ACE_parental_death","ACE_financial","ACE_homelessness","ACE_community_disadvantage",
                  "ACE_community_safety","ACE_peer","ACE_discrimination","ACE_self_illness","ACE_others_illness","ACE_others_death","ACE_natural_disaster","ACE_other_trauma","ACE_unclear", "ACE_original_7", "ACE_non_ACE")
(ACEs_table<-CreateCatTable(vars=ACEs_for_table,data=Scoping_primary, test=F))
print(ACEs_table, showAllLevels = FALSE)


#Funding information reported in the articles

explore_categ(Scoping_primary$Funding_any_info)
explore_categ(Scoping_primary$Funding_number)
682/length(which(Scoping_primary$Funding_any_info=="yes"))
148/length(which(Scoping_primary$Funding_any_info=="yes"))
48/length(which(Scoping_primary$Funding_any_info=="yes"))
14/length(which(Scoping_primary$Funding_any_info=="yes"))

Funding_for_table<-c("Funding_gov", "Funding_academic", "Funding_professional","Funding_private", "Funding_other")
(Funding_table<-CreateCatTable(vars=Funding_for_table,data=Scoping_primary, test=F))
print(Funding_table, showAllLevels = FALSE)

583/length(which(Scoping_primary$Funding_any_info=="yes")) #government funding
154/length(which(Scoping_primary$Funding_any_info=="yes")) #academic funding
38/length(which(Scoping_primary$Funding_any_info=="yes")) #professional society funding
242/length(which(Scoping_primary$Funding_any_info=="yes")) #private funding
71/length(which(Scoping_primary$Funding_any_info=="yes")) #other/unclear

Funding_sources_for_table<-c("Funding_NCI",	"Funding_NHLBI",	"Funding_NIA",	"Funding_NIAAA",	"Funding_NIAID",	"Funding_NIAMS",	"Funding_NIBIB",	"Funding_NICHD",	"Funding_NIDCD",	"Funding_NIDCR",	"Funding_NIDDK",	"Funding_NIDA",	"Funding_NIEHS",	"Funding_NIGMS",	"Funding_NIMH",	"Funding_NIMHD",	"Funding_NINDS",	"Funding_NINR",	"Funding_NLM",	"Funding_Fogarty",	"Funding_NCATS",	"Funding_NCCIH",	"Funding_ORWH",	"Funding_OBSSR","Funding_NIH",	"Funding_ACF",	"Funding_AHRQ",	"Funding_CDC",	"Funding_CMS",	"Funding_HRSA",	"Funding_SAMHSA",	"Funding_VA",	"Funding_NSF",	"Funding_CIHR",	"Funding_NHS",	"Funding_MRC",	"Funding_WHO",	"Funding_Garfield","Funding_Wellcome",	"Funding_WilliamTGrant",	"Funding_RWJ",	"Funding_MacArthur",	"Funding_Jacobs",	"Funding_Casey",	"Funding_NARSAD",	"Funding_Novartis",	"Funding_Janssen",	"Funding_Glaxo",	"Funding_Lilly",	"Funding_Pfizer",	"Funding_Squibb",	"Funding_Astra",	"Funding_Preventive")
Funding_table_by_source<-CreateCatTable(vars=Funding_sources_for_table, data=Scoping_primary, test=F)
print(Funding_table_by_source, showAllLevels = FALSE)

CrossTable(Scoping_primary$Year, Scoping_primary$Funding_ACF, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_AHRQ, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_CDC, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_CMS, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_HRSA, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_NIH, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 
CrossTable(Scoping_primary$Year, Scoping_primary$Funding_SAMHSA, digits=2,chisq=TRUE, prop.chisq=F, prop.r  = T, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS") 

Funding_agency_data<- subset(ScopingR, ScopingR$Funding_agency_org_chart != "" & ScopingR$Funding_year_org_chart !="" & ScopingR$Funding_agency_number_org_chart != "", select=c("Funding_agency_org_chart",  "Funding_year_org_chart","Funding_agency_number_org_chart"))

(Figure4<-ggplot(Funding_agency_data, aes(x = Funding_year_org_chart, y = Funding_agency_number_org_chart)) + 
    stat_summary(fun=identity, geom="point")+stat_summary(fun = mean, geom = "line",aes(group = Funding_agency_org_chart), linetype="solid",size = 1)+
    geom_point(aes(shape = Funding_agency_org_chart), size = 4)+
    scale_shape_manual(values=c(0, 2, 1, 7, 8, 16, 3))+
    labs(x = "Year of publication", y = "Number of articles/year", shape = "Funding agency") +
    scale_x_continuous(breaks=seq(1999,2019,2))+
    scale_y_continuous(breaks=seq(0, 80, by=10))+
    theme(legend.title = element_text (size=18))+
    theme(legend.text=element_text(size=16),legend.position=c(0.15, 0.80))+
    theme(axis.text=element_text(color="black",size=16))+
    theme(axis.title.x = element_text (size=18, vjust = -4))+
    theme(axis.title.y = element_text (size=18, vjust = 5))+
    theme(panel.grid.major = element_line(colour="gray",linetype="dotted", size=0.1),panel.background = element_rect(fill="white") ,panel.border = element_rect(color = "black", fill = NA, size = 1)))+
  theme(plot.background=element_rect(fill="white"),plot.margin = unit(c(1, 1, 1, 1), "cm") )

ggsave("Figure4.jpeg", dpi=800, height=8, width=9.5, units="in")       

explore_categ(Scoping_primary$Funding_conflict_statement)
explore_categ(Scoping_primary$Funding_conflict_declaration)
67/812 # divided by the # of studies in which there was a conflict statement

explore_categ(Scoping_primary$Pharmaceutical_conflict); 27/67


#Country of the samples
explore_categ(Scoping_primary$Country) 
776/1281 #US
69/1281 #UK
64/1281 #Canada
34/1281 #Australia

38/1281 #multinational
12/1281 #unclear

length(unique(Scoping_primary$Country))-6 # of unique entries, excluding multiple and unclear & US, UK,Canada, Australia = number of countries for the rest of the samples
1281-776-69-64-34-38-12 # of articles from other countries
288/1281 # % of articles from other countries


mapped_data <- joinCountryData2Map(Scoping_primary, joinCode = "ISO3",nameJoinColumn = "Unique_country_name",nameCountryColumn = "TableName", verbose = TRUE)
par(mai=c(1,0,0.8,0),xaxs="i",yaxs="i")
(MyMap<-mapCountryData(mapped_data, catMethod = "pretty", numCats = 10, nameColumnToPlot = "No_articles_country", mapTitle="", oceanCol="white", missingCountryCol = "white", addLegend=FALSE, borderCol="black", colourPalette = "white2Black"))
do.call(addMapLegend,c(MyMap,legendLabels="limits",legendWidth=1,legendIntervals="data",legendMar = 4))#adds a legend


#Spin off (is this study affiliated with other studies or surveillance efforts?)

explore_categ(Scoping_primary$Spin_off)
explore_categ(Scoping_primary$Spin_off_or_not_f)
explore_categ(Scoping_primary$Spin_off_name)
explore_categ(Scoping_primary$Spin_off_main_studies)

86*100/1281 #BRFSS studies 
64*100/1281 #NSCH studies
48*100/1281 #Kaiswer studies
21*100/1281 #NESARC studies
20*100/1281 #Florida DJJ studies
1281*.01 # 1% threshold

(86 + 64 + 48 + 21 + 20)/1281 #total number of articles by the top 5 sources

Spin_off_data<- subset(ScopingR, ScopingR$Spin_off_year != "" & ScopingR$Spin_off_category != "" & ScopingR$Spin_off_number !="" & ScopingR$Spin_off_proportion_year , select=c("Spin_off_year","Spin_off_category", "Spin_off_number", "Spin_off_proportion_year"))

(Figure5<-ggplot(Spin_off_data, aes(x = Spin_off_year, y = Spin_off_proportion_year, fill = Spin_off_category, label = Spin_off_number)) +
    geom_bar(stat = "identity") +
    geom_text(size = 5.5, position = position_stack(vjust = 0.5))+
    labs(x = "Year of publication", y = "Proportion of articles/year")+
    scale_x_continuous(breaks=seq(1999,2019,2))+
    scale_y_continuous(breaks=seq(0, 1.0, by=0.1))+
    scale_fill_manual(values = c("yes" = "gray55", "no" = "gray85"))+
    guides(fill=guide_legend(title="Part of a larger study?"))+
    theme(legend.title = element_text (size=20))+
    theme(legend.text=element_text(size=18)))+
  theme(axis.text=element_text(color="black",size=16))+
  theme(axis.title.x = element_text (color = "black",  size = 20,vjust=-4))+
  theme(axis.title.y = element_text(color = "black", size=20, vjust=5)) +
  theme(panel.grid.major = element_line(colour="white",linetype="dashed", size=0.1),panel.background = element_rect(fill="white") ,panel.border = element_rect(color = "black", fill = NA))+
  theme(plot.background=element_rect(fill="white"),plot.margin = unit(c(1, 2, 1, 2), "cm"))

ggsave("Figure5.jpeg", dpi=800, height=10, width=14, units="in")   


#How many articles were eligible for data analyses for the WHO study?
(length(which(Scoping_primary$WHO_eligible=="yes")))

#How many articles were coded as unclear/inconsistent?
(length(which(Scoping_primary$WHO_unclear=="yes")))
(length(which(Scoping_primary$WHO_unclear=="yes")))/(length(which(Scoping_primary$WHO_eligible=="yes"))) # proportion of articles coded as unclear/inconsistent

#How many articles were coded, excluding unclear articles?
(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

#What were the variable roles assigned to each category in the WHO framework? 

explore_categ(Scoping_primary$WHO_governance)
explore_categ(Scoping_primary$WHO_culture)
explore_categ(Scoping_primary$WHO_social_class)
explore_categ(Scoping_primary$WHO_wealth)
explore_categ(Scoping_primary$WHO_income)
explore_categ(Scoping_primary$WHO_education)
explore_categ(Scoping_primary$WHO_occupation)
explore_categ(Scoping_primary$WHO_gender)
explore_categ(Scoping_primary$WHO_race_ethnicity)
explore_categ(Scoping_primary$WHO_geographic)
explore_categ(Scoping_primary$WHO_other_social_determinant)
explore_categ(Scoping_primary$WHO_clear_composite)
explore_categ(Scoping_primary$WHO_unclear_composite)
explore_categ(Scoping_primary$WHO_welfare_services)
explore_categ(Scoping_primary$WHO_housing)
explore_categ(Scoping_primary$WHO_hunger)
explore_categ(Scoping_primary$WHO_neighborhood)
explore_categ(Scoping_primary$WHO_other_material)
explore_categ(Scoping_primary$WHO_household_composition)
explore_categ(Scoping_primary$WHO_parental_psychosocial)
explore_categ(Scoping_primary$WHO_ACEs)
explore_categ(Scoping_primary$WHO_psychosocial)
explore_categ(Scoping_primary$WHO_other_psychosocial)
explore_categ(Scoping_primary$WHO_health_risk_promotive)
explore_categ(Scoping_primary$WHO_genetic)
explore_categ(Scoping_primary$WHO_biomarkers)
explore_categ(Scoping_primary$WHO_health_system)
explore_categ(Scoping_primary$WHO_child_welfare_system)
explore_categ(Scoping_primary$WHO_bio_health)
explore_categ(Scoping_primary$WHO_nonbio_health)
explore_categ(Scoping_primary$WHO_mortality)
explore_categ(Scoping_primary$WHO_justice_system)
explore_categ(Scoping_primary$WHO_parenthood)
explore_categ(Scoping_primary$WHO_unclear)

WHO_for_table<-c(	"WHO_governance","WHO_culture","WHO_social_class","WHO_wealth", "WHO_income", "WHO_education", "WHO_occupation", "WHO_gender","WHO_race_ethnicity", "WHO_geographic", "WHO_other_social_determinant", "WHO_clear_composite", "WHO_unclear_composite", "WHO_welfare_services", "WHO_housing", "WHO_hunger", "WHO_neighborhood", "WHO_other_material", "WHO_household_composition", "WHO_parental_psychosocial", "WHO_ACEs", "WHO_psychosocial", "WHO_other_psychosocial", "WHO_health_risk_promotive", "WHO_genetic", "WHO_biomarkers", "WHO_health_system", "WHO_child_welfare_system", "WHO_bio_health", "WHO_nonbio_health", "WHO_mortality", "WHO_justice_system", "WHO_parenthood", "WHO_unclear")
WHO_table<-CreateCatTable(vars=WHO_for_table,data=Scoping_primary, test=F)
print(WHO_table, showAllLevels = FALSE)

#How many articles had at least something coded as bidirectional?
explore_categ(Scoping_primary$WHO_bidirectional)
10/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))*100


#Which variables were coded as population, predictor, moderator, and grouping variables?

WHO_for_population_table<-c("WHO_governance_population", "WHO_culture_population","WHO_social_class_population","WHO_wealth_population","WHO_income_population","WHO_education_population",	"WHO_occupation_population","WHO_gender_population","WHO_race_ethnicity_population", "WHO_geographic_population","WHO_other_social_determinant_population","WHO_clear_composite_population","WHO_unclear_composite_population","WHO_welfare_services_population","WHO_housing_population","WHO_hunger_population","WHO_neighborhood_population","WHO_other_material_population","WHO_household_composition_population",	"WHO_parental_psychosocial_population","WHO_parenthood_population","WHO_ACEs_population",	"WHO_psychosocial_population","WHO_other_psychosocial_population","WHO_health_risk_promotive_population",	"WHO_genetic_population","WHO_biomarkers_population","WHO_health_system_population",	"WHO_child_welfare_system_population","WHO_bio_health_population","WHO_nonbio_health_population","WHO_mortality_population","WHO_justice_system_population")
WHO_population_table<-CreateCatTable(vars=WHO_for_population_table,data=Scoping_primary, test=F)
print(WHO_population_table, showAllLevels = FALSE)
(84/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear"))))*100

WHO_for_predictor_table<-c("WHO_governance_predictor", "WHO_culture_predictor","WHO_social_class_predictor","WHO_wealth_predictor","WHO_income_predictor","WHO_education_predictor",	"WHO_occupation_predictor","WHO_gender_predictor","WHO_race_ethnicity_predictor", "WHO_geographic_predictor","WHO_other_social_determinant_predictor","WHO_clear_composite_predictor","WHO_unclear_composite_predictor","WHO_welfare_services_predictor","WHO_housing_predictor","WHO_hunger_predictor","WHO_neighborhood_predictor","WHO_other_material_predictor","WHO_household_composition_predictor",	"WHO_parental_psychosocial_predictor","WHO_parenthood_predictor","WHO_ACEs_predictor",	"WHO_psychosocial_predictor","WHO_other_psychosocial_predictor","WHO_health_risk_promotive_predictor",	"WHO_genetic_predictor","WHO_biomarkers_predictor","WHO_health_system_predictor",	"WHO_child_welfare_system_predictor","WHO_bio_health_predictor","WHO_nonbio_health_predictor","WHO_mortality_predictor","WHO_justice_system_predictor")
WHO_predictor_table<-CreateCatTable(vars=WHO_for_predictor_table,data=Scoping_primary, test=F)
print(WHO_predictor_table, showAllLevels = FALSE)
62/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

WHO_for_grouping_table<-c("WHO_governance_grouping", "WHO_culture_grouping","WHO_social_class_grouping","WHO_wealth_grouping","WHO_income_grouping","WHO_education_grouping",	"WHO_occupation_grouping","WHO_gender_grouping","WHO_race_ethnicity_grouping", "WHO_geographic_grouping","WHO_other_social_determinant_grouping","WHO_clear_composite_grouping","WHO_unclear_composite_grouping","WHO_welfare_services_grouping","WHO_housing_grouping","WHO_hunger_grouping","WHO_neighborhood_grouping","WHO_other_material_grouping","WHO_household_composition_grouping",	"WHO_parental_psychosocial_grouping","WHO_parenthood_grouping","WHO_ACEs_grouping",	"WHO_psychosocial_grouping","WHO_other_psychosocial_grouping","WHO_health_risk_promotive_grouping",	"WHO_genetic_grouping","WHO_biomarkers_grouping","WHO_health_system_grouping",	"WHO_child_welfare_system_grouping","WHO_bio_health_grouping","WHO_nonbio_health_grouping","WHO_mortality_grouping","WHO_justice_system_grouping")
WHO_grouping_table<-CreateCatTable(vars=WHO_for_grouping_table,data=Scoping_primary, test=F)
print(WHO_grouping_table, showAllLevels = FALSE)
62/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

WHO_for_moderator_table<-c("WHO_governance_moderator", "WHO_culture_moderator","WHO_social_class_moderator","WHO_wealth_moderator","WHO_income_moderator","WHO_education_moderator",	"WHO_occupation_moderator","WHO_gender_moderator","WHO_race_ethnicity_moderator", "WHO_geographic_moderator","WHO_other_social_determinant_moderator","WHO_clear_composite_moderator","WHO_unclear_composite_moderator","WHO_welfare_services_moderator","WHO_housing_moderator","WHO_hunger_moderator","WHO_neighborhood_moderator","WHO_other_material_moderator","WHO_household_composition_moderator",	"WHO_parental_psychosocial_moderator","WHO_parenthood_moderator","WHO_ACEs_moderator",	"WHO_psychosocial_moderator","WHO_other_psychosocial_moderator","WHO_health_risk_promotive_moderator",	"WHO_genetic_moderator","WHO_biomarkers_moderator","WHO_health_system_moderator",	"WHO_child_welfare_system_moderator","WHO_bio_health_moderator","WHO_nonbio_health_moderator","WHO_mortality_moderator","WHO_justice_system_moderator")
WHO_moderator_table<-CreateCatTable(vars=WHO_for_moderator_table,data=Scoping_primary, test=F)
print(WHO_moderator_table, showAllLevels = FALSE)
87/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))


#What was coded as mediator?
WHO_for_mediator_table<-c("WHO_governance_mediator", "WHO_culture_mediator","WHO_social_class_mediator","WHO_wealth_mediator","WHO_income_mediator","WHO_education_mediator",	"WHO_occupation_mediator","WHO_gender_mediator","WHO_race_ethnicity_mediator", "WHO_geographic_mediator","WHO_other_social_determinant_mediator","WHO_clear_composite_mediator","WHO_unclear_composite_mediator","WHO_welfare_services_mediator","WHO_housing_mediator","WHO_hunger_mediator","WHO_neighborhood_mediator","WHO_other_material_mediator","WHO_household_composition_mediator",	"WHO_parental_psychosocial_mediator","WHO_parenthood_mediator","WHO_ACEs_mediator",	"WHO_psychosocial_mediator","WHO_other_psychosocial_mediator","WHO_health_risk_promotive_mediator",	"WHO_genetic_mediator","WHO_biomarkers_mediator","WHO_health_system_mediator",	"WHO_child_welfare_system_mediator","WHO_bio_health_mediator","WHO_nonbio_health_mediator","WHO_mortality_mediator","WHO_justice_system_mediator")
WHO_mediator_table<-CreateCatTable(vars=WHO_for_mediator_table,data=Scoping_primary, test=F)
print(WHO_mediator_table, showAllLevels = FALSE)
78/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))


#How was health coded across all articles?
WHO_for_health_table<-c("WHO_health_population",	"WHO_health_predictor",	"WHO_health_grouping",	"WHO_health_moderator",	"WHO_health_mediator",	"WHO_health_outcome",	"WHO_health_bidirectional")
WHO_health_table<-CreateCatTable(vars=WHO_for_health_table,data=Scoping_primary, test=F)
print(WHO_health_table, showAllLevels = FALSE)
(51/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear"))))*100


#Was there any reference to social determinants of health inequities in the research goals related to ACEs?
explore_categ(Scoping_primary$WHO_anything_social)
655/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))*100

#Aggregate results within the social & intermediary levels, i.e., how many articles within the "social determinants of health inequities" & "determinants of health" & "health" (i.e., bio, non-bio, mortality) levels had something coded as population, etc. in the research goals?
explore_categ(Scoping_primary$WHO_social_population)
explore_categ(Scoping_primary$WHO_social_predictor)
explore_categ(Scoping_primary$WHO_social_grouping)
explore_categ(Scoping_primary$WHO_social_moderator)
explore_categ(Scoping_primary$WHO_social_mediator)
explore_categ(Scoping_primary$WHO_social_health)
explore_categ(Scoping_primary$WHO_social_bidirectional)

explore_categ(Scoping_primary$WHO_intermediary_population)
explore_categ(Scoping_primary$WHO_intermediary_predictor)
explore_categ(Scoping_primary$WHO_intermediary_grouping)
explore_categ(Scoping_primary$WHO_intermediary_moderator)
explore_categ(Scoping_primary$WHO_intermediary_mediator)
explore_categ(Scoping_primary$WHO_intermediary_health)
explore_categ(Scoping_primary$WHO_intermediary_bidirectional)

explore_categ(Scoping_primary$WHO_health_population)
explore_categ(Scoping_primary$WHO_health_predictor)
explore_categ(Scoping_primary$WHO_health_grouping)
explore_categ(Scoping_primary$WHO_health_moderator)
explore_categ(Scoping_primary$WHO_health_mediator)
explore_categ(Scoping_primary$WHO_health_outcome)
explore_categ(Scoping_primary$WHO_health_bidirectional)

4/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))*100

#What was included in the SES composites?
explore_categ(Scoping_primary$WHO_SES_composite)
32/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))


#Within the intermediary determinants level, which domains were included in how many research goals?
explore_categ(Scoping_primary$WHO_intermediary_material)
explore_categ(Scoping_primary$WHO_intermediary_psychosocial)
explore_categ(Scoping_primary$WHO_intermediary_biobehavioral)
explore_categ(Scoping_primary$WHO_intermediary_systems)

#What were the roles assigned to the health system vs. child welfare vs. justice systems in the research goals?
explore_categ(Scoping_primary$WHO_health_system)
explore_categ(Scoping_primary$WHO_child_welfare_system)
explore_categ(Scoping_primary$WHO_justice_system)
(24/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear"))))*100


#Articles with biological vs. nonbiological health outcomes
explore_categ(Scoping_primary$WHO_any_nonbio_outcome)
explore_categ(Scoping_primary$WHO_any_bio_outcome)
(310/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear"))))*100


#How many articles included only predictor-health(Story of Decline)?

explore_categ(Scoping_primary$WHO_predictor_outcome_only)
592/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

#How many articles included only predictor-moderator/mediator-health (Story of Helplessness & Control)?
explore_categ(Scoping_primary$WHO_predictor_moderator_mediator_outcome)
471/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

#How many articles fit the WHO Story? 
explore_categ(Scoping_primary$WHO_consistent_w_WHO_story)
307/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

explore_categ(Scoping_primary$WHO_consistent_w_WHO_story_excluding_population)
117/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

#How many articles fit the Intermediary Determinants Story? 
explore_categ(Scoping_primary$WHO_consistent_w_intermediary_determinant_story)
784/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))

#How many articles included intermediary determinants as predictor & social determinants as health (the Intergenerational Transmission Story)?
explore_categ(Scoping_primary$WHO_consistent_w_intergenerational_transmission_story)
47/(length(which(Scoping_primary$WHO_eligible_and_clear=="eligible and clear")))


#Protective Factors
explore_categ(Scoping_primary$Protective_any)

Protective_for_table<-c("Protective_individual", "Protective_family_psychosocial", "Protective_family_financial","Protective_neighborhood", "Protective_societal", "Protective_unclear")
(Protective_table<-CreateCatTable(vars=Protective_for_table,data=Scoping_primary, test=F))
print(Protective_table, showAllLevels = FALSE)

69/142 #individual protective factors/# of studies that had a protective factor
58/142 #family-psychosocial/# of studies that had a protective factor
6/142 #family-financial/# of studies that had a protective factor
80/142 #neighborhood/# of studies that had a protective factor
6/142 #societal/# of studies that had a protective factor
86/142 #societal + community/neighborhood/# of studies that had a protective factor


CrossTable(Scoping_primary$Year, Scoping_primary$Protective_any, digits=0,chisq=TRUE, prop.chisq=F, prop.r  = F, prop.t = F, fisher=F, expected = F, sresid=FALSE, format="SPSS")  
(27+33+48)/142


#Recommendations

Recommendations_for_table<-c(	"Recommendations_primary_primordial","Recommendations_secondary_tertiary","Recommendations_both","Recommendations_vague_none")
Recommendations_table<-CreateCatTable(vars=Recommendations_for_table,data=Scoping_primary, test=F)
print(Recommendations_table, showAllLevels = FALSE)

(12+208)/1282
(547+208)/1282

(10/700)*100
(267/700)*100
(107/700)*100
(248/700)*100

Recommendations_table_by_year<-CreateCatTable(vars=Recommendations_for_table,strata = c("Year"), data=Scoping_primary, test=F)
print(Recommendations_table_by_year, strata = c("Year"), showAllLevels = FALSE)

Recommendations_data<- subset(ScopingR, ScopingR$Recommendation_year !="NA" & ScopingR$Recommendation_year_type !="NA" & ScopingR$Recommendation_per_year!= "NA", ScopingR$Recommendation_proportion_year!="NA", ScopingR$Recommendation_year_number !="NA", select=c("Recommendation_year", "Recommendation_year_type", "Recommendation_per_year", "Recommendation_proportion_year", "Recommendation_year_number"))

(Figure7<-ggplot(Recommendations_data, aes(x = Recommendation_year, y = Recommendation_proportion_year, fill = Recommendation_year_type, label = Recommendation_per_year)) +
    geom_bar(stat = "identity") +
    geom_text(size = 5.5, position = position_stack(vjust = 0.5))+
    labs(x = "Year of publication", y = "Proportion of articles")+
    scale_x_continuous(breaks=seq(1999,2019,2))+
    scale_y_continuous(breaks=seq(0, 1.0, by=0.1))+
    scale_fill_manual(values = c("Both" = "gray50", "Primordial or Primary" = "gray65","Secondary or Tertiary"="gray80","Vague or None"="gray95"))+
    guides(fill=guide_legend(title="Recommendation Type"))+
    theme(legend.title = element_text (size=20))+
    theme(legend.text=element_text(size=18)))+
  theme(axis.text=element_text(color="black",size=16))+
  theme(axis.title.x = element_text (color = "black",  size = 20,vjust=-4))+
  theme(axis.title.y = element_text(color = "black", size=20, vjust=5)) +
  theme(panel.grid.major = element_line(colour="white",linetype="dashed", size=0.1),panel.background = element_rect(fill="white") ,panel.border = element_rect(color = "black", fill = NA))+
  theme(plot.background=element_rect(fill="white"),plot.margin = unit(c(1, 2, 1, 2), "cm"))

ggsave("Figure7.jpeg", dpi=800, height=10, width=14, units="in")