---
title: "Figure 2 source code eLife submission"
output: html_document
---

```{r message = FALSE, warning = FALSE}
library(easyGgplot2)
#ggplot2.multiplot(..., plotlist=NULL, cols=2)
library(easyGgplot2)
library(lme4)
library(merTools)
combinedData <- read.csv("~/Desktop/Temperature study/combined_data.csv")

combinedData$time_HR_imputed <- combinedData$time_HR
combinedData$time_HR_imputed[is.na(combinedData$time_HR_imputed)] <- 12

dataVeterans  <- read.csv('~/Desktop/Temperature study/veterans_processed.csv')
dataNHANES    <- read.csv('~/Desktop/Temperature study/NHANES_processed.csv')
ST_data 	  <- read.csv('~/Desktop/Temperature study/STRIDE_processed.csv')

dataVeterans$race <- factor(as.character(dataVeterans$race),levels = c("white", "black"))
combinedData$race [combinedData$race=="asian"]<- "other"
combinedData$race <- factor(combinedData$race,levels = c("white", "black","other" ))
dataNHANES$race[dataNHANES$race == "unknown"] <- NA
dataNHANES$race <- factor(dataNHANES$race,levels = c("white", "black","other" ))
ST_data$race [ST_data$race=="Asian"]<- "Other"
ST_data$race <- factor(ST_data$race,levels = c("White", "Black","Other"))
ST_data$exam_month<-factor (ST_data$exam_month)

dataVeterans$state_type <- factor(dataVeterans$state_type, levels = c("very cold", "cold","moderately_cold","moderate","warm", "hot"))

dataNHANES$exammonth <-as.factor(dataNHANES$exammonth)
dataVeterans$exammonth <-as.factor(dataVeterans$exammonth)


state_temps <- read.csv("~/Desktop/Temperature study/state_temps.csv")
state_codes <- read.csv("~/Desktop/Temperature study/state_codes.csv")
state_codes$Code <- factor(as.character(state_codes$Code), levels=levels(dataVeterans$state))
state_temps$state_code <- state_codes$Code[match(as.character(state_temps$state),as.character(state_codes$State))]

dataVeterans$ambient_temp <- (state_temps$value[match(paste(dataVeterans$state,(dataVeterans$examyear * 100 + as.numeric(dataVeterans$exammonth))), paste(state_temps$state_code,state_temps$date))] - 32) / 1.8

ST_data$state <- "CA"
  
ST_data$ambient_temp <- (state_temps$value[match(paste(ST_data$state,(ST_data$exam_year * 100 + as.numeric(ST_data$exam_month))), paste(state_temps$state_code,state_temps$date))] - 32) / 1.8


```

Figure 2: Temperature trends within birth cohorts of the UAVCW, 1860-1940 (black and white men). 
Smoothed unadjusted data (local regression) for temperature measurement trends within birth cohorts.  The different colors represent different birth cohorts (green: 1820s, blue: 1830s, orange: 1840s)

```{r message = FALSE, warning = FALSE}
# Minimum ages vary across the three cohorts due to the bias introduced by the sampling (it is impossible to be younger than 30, be born in 1820s, and be a veteran of the Civil War). The values 52, 45, and 38 represent the cutoffs so that 99% of observations in a respective cohort are older than this cutoff.

Veterans_birth_cohorts_C  <-ggplot(data=dataVeterans[ ((dataVeterans$age>52 & dataVeterans$age<80 & dataVeterans$birth_cohort== "1820s") | (dataVeterans$age>45 & dataVeterans$age<80 & dataVeterans$birth_cohort== "1830s") | (dataVeterans$age>38 & dataVeterans$age<80 & dataVeterans$birth_cohort== "1840s"))  & dataVeterans$race == "white", ], aes(age, temp_C, group=birth_cohort, colour= birth_cohort))   + geom_smooth()

Veterans_birth_cohorts_C + theme(legend.position="bottom", axis.text=element_text(size=14), axis.title=element_text(size=14), legend.text=element_text(size=14)) + labs(x= "Age, years", y = "Temperature,°C", colour = "Birth decade" ) + scale_colour_manual(values = plot_colors[c(2,3,1 )]) #+ylim(36.7,37)
```

```{r message = FALSE, warning = FALSE} 
lmBirthcohorts <-lm(temp_C ~ age + weightKG + heightCM + race + birth_cohort, data
 =dataVeterans [dataVeterans$birth_cohort %in%c ( "1820s", "1830s",  "1840s"), ])
 summary (lmBirthcohorts)
```  

And we present predicted temperature for men at age 30, weight 70 kg, and hieght 170 cm for each of these decades in a table

```{r message = FALSE, warning = FALSE}
 predict_dataframe <- data.frame(age=30, weightKG=70, heightCM = 170, race= "white", birth_cohort="1820s")
 round(predict (lmBirthcohorts , predict_dataframe , interval="confidence"), 2)
 
  
predict_dataframe <- data.frame(age=30, weightKG=70, heightCM = 170,race= "white",birth_cohort="1830s")
 round(predict (lmBirthcohorts , predict_dataframe , interval="confidence"),2)
```

