---
title: "Figure supplement 3"
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 1-figure supplement 3. Model of mean body temperature in black and white ethnicity groups in three different time periods (cohorts): 1860-1940 (UAVCW), 1971-1975 (NHANES 1), and 2007-2017 (STRIDE).
A. Coefficients (and standard errors) from multivariate linear regression models for each cohort including age, body weight, height and time of day. Separate regression models were fitted for men and women, in black and white ethnicity groups.  NHANES I includes population weights. *, **, *** indicates significance at the 90%, 95%, and 99% level, respectively. B. Expected body temperature (and associated 95% confidence interval) for 30-year old men and women with body weight 70 kg and height 170 cm in each time period/cohort. These values derive from the regression models presented in A.


```{r message = FALSE, warning = FALSE}
lmVeterans_whites <-lm (temp_C~ age + weightKG + heightCM , data= dataVeterans [dataVeterans$race == "white",])
summary(lmVeterans_whites)
lmVeterans_whites_predict =data.frame(age=30, weightKG=70, heightCM = 170, race="white")

round(predict(lmVeterans_whites, lmVeterans_whites_predict, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}
lmVeterans_blacks <-lm (temp_C~ age + weightKG + heightCM , data= dataVeterans [dataVeterans$race == "black",])
summary(lmVeterans_blacks)

lmVeterans_b_predict =data.frame(age=30, weightKG=70, heightCM = 170, race="black")
round(predict(lmVeterans_blacks, lmVeterans_b_predict, interval="confidence"),2)

```

```{r message = FALSE, warning = FALSE}
lmNHANES_W_HR <- lm(temp_C ~ age +  weight_KG + height_CM + time_HR, data=dataNHANES[dataNHANES$sex == "male" & dataNHANES$race == "white", ], weights=dataNHANES$sample_weights[dataNHANES$sex=="male" & dataNHANES$race == "white"])
summary(lmNHANES_W_HR)

nobs(lmNHANES_W_HR)
NHANES_predict=data.frame (age=30, weight_KG=70, height_CM = 170, time_HR=12, race="white", sex ="male")
round(predict(lmNHANES_W_HR, NHANES_predict, interval="confidence"),2)
```


```{r message = FALSE, warning = FALSE}
lmNHANES_W_black_HR <- lm(temp_C ~ age + weight_KG + height_CM + time_HR, data=dataNHANES[dataNHANES$sex == "male" & dataNHANES$race == "black", ], weights=dataNHANES$sample_weights[dataNHANES$sex=="male" & dataNHANES$race == "black"])
summary(lmNHANES_W_black_HR)
nobs(lmNHANES_W_black_HR)
        
NHANES_predict_black=data.frame (age=30, weight_KG=70, height_CM = 170, race="black", sex ="male", time_HR= 12)
round(predict(lmNHANES_W_black_HR, NHANES_predict_black, interval="confidence"),2)

```

```{r message = FALSE, warning = FALSE}


lmNHANES_W_F_HR <- lm(temp_C ~ age +  weight_KG + height_CM + time_HR, data=dataNHANES[dataNHANES$sex == "female" & dataNHANES$race == "white", ], weights=dataNHANES$sample_weights[dataNHANES$sex=="female" & dataNHANES$race == "white"])
summary(lmNHANES_W_F_HR)
nobs(lmNHANES_W_F_HR)
NHANES_predict=data.frame (age=30, weight_KG=70, height_CM = 170, time_HR=12, race="white", sex ="female")
round(predict(lmNHANES_W_F_HR , NHANES_predict, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}


lmNHANES_W_F_HR_black <- lm(temp_C ~ age +  weight_KG + height_CM + time_HR, data=dataNHANES[dataNHANES$sex == "female" & dataNHANES$race == "black", ], weights=dataNHANES$sample_weights[dataNHANES$sex=="female" & dataNHANES$race == "black"])
summary(lmNHANES_W_F_HR_black)
nobs(lmNHANES_W_F_HR_black)
NHANES_predict_b_f=data.frame (age=30, weight_KG=70, height_CM = 170, time_HR=12, race="white", sex ="female")
round(predict(lmNHANES_W_F_HR_black , NHANES_predict_b_f, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}

 lmSTRIDE_HR <- lm(temp_C ~ age_years +  weight_KG + height_CM + time_HR, data=ST_data[ST_data$GENDER == "Male" & ST_data$race == "White", ])
summary(lmSTRIDE_HR)
nobs(lmSTRIDE_HR)
STRIDE_predict=data.frame(age_years=30, weight_KG=70, height_CM = 170, time_HR=12, race="White", GENDER = "Male")

round(predict(lmSTRIDE_HR, STRIDE_predict, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}

 lmSTRIDE_HR_B <- lm(temp_C ~ age_years +  weight_KG + height_CM + time_HR, data=ST_data[ST_data$GENDER == "Male" & ST_data$race == "Black", ])
summary(lmSTRIDE_HR_B)
nobs(lmSTRIDE_HR_B)
STRIDE_predict_B=data.frame(age_years=30, weight_KG=70, height_CM = 170, time_HR=12, race="Black", GENDER = "Male")

round(predict(lmSTRIDE_HR_B, STRIDE_predict_B, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}

 lmSTRIDE_HR_B_F <- lm(temp_C ~ age_years +  weight_KG + height_CM + time_HR, data=ST_data[ST_data$GENDER == "Female" & ST_data$race == "Black", ])
summary(lmSTRIDE_HR_B_F)
nobs(lmSTRIDE_HR_B_F)
STRIDE_predict_B_F=data.frame(age_years=30, weight_KG=70, height_CM = 170, time_HR=12, race="Black", GENDER = "Female")

round(predict(lmSTRIDE_HR_B_F, STRIDE_predict_B_F, interval="confidence"),2)
```

```{r message = FALSE, warning = FALSE}

 lmSTRIDE_F_HR <- lm(temp_C ~ age_years + weight_KG + height_CM + time_HR ,data=ST_data[ST_data$GENDER == "Female" & ST_data$race == "White", ])
summary(lmSTRIDE_F_HR)

STRIDE_predict=data.frame(age_years=30, weight_KG=70, height_CM = 170, time_HR=12, race="White", GENDER = "Female")

round(predict(lmSTRIDE_F_HR, STRIDE_predict, interval="confidence"),2)
```


