Loading dependent variables
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
setwd("C:\\Users\\anami\\OneDrive\\Documents\\Health Disparities\\Project PA")
data1<- read.csv("C:\\Users\\anami\\OneDrive\\Documents\\Health Disparities\\Project PA\\data PLACES.csv")
Independent variable: Walkability Index
dataA <- read.csv("C:\\Users\\anami\\OneDrive\\Documents\\Health Disparities\\Project PA\\WalkabilityIndex_Tract_2019.csv" )
library(dplyr)
dataB<-dataA%>%
filter(StAbbr=="TX")
dataB
Calculating walkability index for Texas census tracts
library(ggplot2)
dataB$rank1 <- as.numeric(cut_number(dataB$D2A_EPHHM,20))
dataB$rank2 <- as.numeric(cut_number(dataB$D2B_E8MIXA,20))
dataB$rank3 <- as.numeric(cut_number(dataB$D3B,20))
library(dplyr)
dataC<-subset(dataB, select=c(TractID2019,Pop2018,rank1,rank2,rank3,D4A,NatWalkInd)) %>%
mutate(d4A =ifelse(D4A ==-99999.00,1, D4A))
dataC
library(ggplot2)
dataB$rank1 <- as.numeric(cut_number(dataB$D2A_EPHHM,20))
dataB$rank2 <- as.numeric(cut_number(dataB$D2B_E8MIXA,20))
dataB$rank3 <- as.numeric(cut_number(dataB$D3B,20))
library(dplyr)
dataC<-subset(dataB, select=c(TractID2019,Pop2018,rank1,rank2,rank3,D4A,NatWalkInd)) %>%
mutate(d4A =ifelse(D4A ==-99999.00,1, D4A))
dataC
dataC$rank4 <- as.numeric( cut(dataC$d4A,20, labels=F))
library(dplyr)
dataD<-dataC %>%
mutate(WI=(rank1/6)+(rank2/6)+(rank3/3)+(rank4/3))
names(dataD)
## [1] "TractID2019" "Pop2018" "rank1" "rank2" "rank3"
## [6] "D4A" "NatWalkInd" "d4A" "rank4" "WI"
library(dplyr)
dataE<-subset(dataD, select=c(TractID2019,WI)) %>%
group_by(TractID2019)
names(dataE)
## [1] "TractID2019" "WI"
Control Variables
file1<- read.csv("C:\\Users\\anami\\OneDrive\\Documents\\Health Disparities\\Project PA\\ACS2019.csv")
file2<-subset(file1, select=c(geoid,totpop,fpop,mage,hispan,nh_white,nh_black,nh_asian))
names(file2)
## [1] "geoid" "totpop" "fpop" "mage" "hispan" "nh_white" "nh_black"
## [8] "nh_asian"
nrow(file2)
## [1] 5265
file3<- read.csv("C:\\Users\\anami\\OneDrive\\Documents\\Health Disparities\\Project PA\\SVI2.csv")
file4<-subset(file3, select=c(FIPS,RPL_THEMES))
names(file4)
## [1] "FIPS" "RPL_THEMES"
nrow(file4)
## [1] 5254
summary(file4)
## FIPS RPL_THEMES
## Min. :4.800e+10 Min. :-999.0000
## 1st Qu.:4.811e+10 1st Qu.: 0.2435
## Median :4.820e+10 Median : 0.4957
## Mean :4.823e+10 Mean : -8.0606
## 3rd Qu.:4.836e+10 3rd Qu.: 0.7478
## Max. :4.851e+10 Max. : 1.0000
nrow(file4)
## [1] 5254
library(dplyr)
file5<-file4 %>%
mutate(SVI =ifelse(RPL_THEMES ==-999.000,NA, RPL_THEMES))
file5
file6<-subset(file5, select=c(FIPS,SVI))
names(file6)
## [1] "FIPS" "SVI"
nrow(file6)
## [1] 5254
file7<-file6%>%
na.omit(file6)
nrow(file7)
summary(file7)
## FIPS SVI
## Min. :4.800e+10 Min. :0.00
## 1st Qu.:4.811e+10 1st Qu.:0.25
## Median :4.820e+10 Median :0.50
## Mean :4.823e+10 Mean :0.50
## 3rd Qu.:4.836e+10 3rd Qu.:0.75
## Max. :4.851e+10 Max. :1.00
nrow(file7)
## [1] 5209
Combining datasets:
combineA <-data2 %>%
left_join(dataE, by = c("TractFIPS" = "TractID2019")) %>%
dplyr::select(TractFIPS,WI,OBESITY_CrudePrev,LPA_CrudePrev)
names(combineA)
## [1] "TractFIPS" "WI" "OBESITY_CrudePrev"
## [4] "LPA_CrudePrev"
combineB <-combineA %>%
left_join(file2, by = c("TractFIPS" = "geoid")) %>%
dplyr::select(TractFIPS,WI,OBESITY_CrudePrev,LPA_CrudePrev,totpop,fpop,mage,hispan,nh_white,nh_black,nh_asian)
names(combineB)
## [1] "TractFIPS" "WI" "OBESITY_CrudePrev"
## [4] "LPA_CrudePrev" "totpop" "fpop"
## [7] "mage" "hispan" "nh_white"
## [10] "nh_black" "nh_asian"
combine1 <-combineB %>%
left_join(file7, by = c("TractFIPS" = "FIPS")) %>%
dplyr::select(TractFIPS,WI,OBESITY_CrudePrev, LPA_CrudePrev,totpop,fpop,mage,hispan,nh_white,nh_black,nh_asian,SVI)
names(combine1)
## [1] "TractFIPS" "WI" "OBESITY_CrudePrev"
## [4] "LPA_CrudePrev" "totpop" "fpop"
## [7] "mage" "hispan" "nh_white"
## [10] "nh_black" "nh_asian" "SVI"
nrow(combine1)
## [1] 5222
summary(combine1)
## TractFIPS WI OBESITY_CrudePrev LPA_CrudePrev
## Min. :4.800e+10 Min. : 1.000 Min. :18.30 Min. :10.10
## 1st Qu.:4.811e+10 1st Qu.: 6.000 1st Qu.:32.73 1st Qu.:20.80
## Median :4.820e+10 Median : 8.333 Median :36.65 Median :26.50
## Mean :4.823e+10 Mean : 8.358 Mean :36.83 Mean :27.37
## 3rd Qu.:4.836e+10 3rd Qu.:10.667 3rd Qu.:40.90 3rd Qu.:33.30
## Max. :4.851e+10 Max. :19.000 Max. :55.30 Max. :56.30
##
## totpop fpop mage hispan
## Min. : 16 Min. : 0 Length:5222 Length:5222
## 1st Qu.: 3384 1st Qu.: 1682 Class :character Class :character
## Median : 4812 Median : 2414 Mode :character Mode :character
## Mean : 5412 Mean : 2724
## 3rd Qu.: 6640 3rd Qu.: 3350
## Max. :72041 Max. :37626
##
## nh_white nh_black nh_asian SVI
## Length:5222 Length:5222 Length:5222 Min. :0.0002
## Class :character Class :character Class :character 1st Qu.:0.2501
## Mode :character Mode :character Mode :character Median :0.5002
## Mean :0.5001
## 3rd Qu.:0.7501
## Max. :1.0000
## NA's :15
Variables description
library(ggplot2)
combine1$group <- as.numeric(cut_number(combine1$WI, 4))
combine1$group<-factor(combine1$group, levels = c(1, 2, 3,4),
labels = c("Q1", "Q2", "Q3","Q4"))
table(combine1$group)
##
## Q1 Q2 Q3 Q4
## 1332 1346 1246 1298
combine1a <- combine1 %>%
rename(
LPA = LPA_CrudePrev,
Obesity = OBESITY_CrudePrev)
library(dplyr)
final1<-subset(combine1, select=c(WI,totpop,fpop,mage,hispan,nh_white,nh_black,nh_asian,group))%>%
mutate(percentfemale=(fpop*100)/totpop)
names(final1)
## [1] "WI" "totpop" "fpop" "mage"
## [5] "hispan" "nh_white" "nh_black" "nh_asian"
## [9] "group" "percentfemale"
total_sum <- sum(final1$totpop)
print(total_sum)
## [1] 28260633
library(dplyr)
final1<-subset(combine1, select=c(WI,totpop,fpop,mage,hispan,nh_white,nh_black,nh_asian,SVI,group))%>%
mutate(percentfemale=(fpop*100)/totpop)
names(final1)
## [1] "WI" "totpop" "fpop" "mage"
## [5] "hispan" "nh_white" "nh_black" "nh_asian"
## [9] "SVI" "group" "percentfemale"
final1$mage <- as.numeric(final1$mage)
final1$hispan <- as.numeric(final1$hispan)
final1$nh_white <- as.numeric(final1$nh_white)
final1$nh_black <- as.numeric(final1$nh_black)
final1$nh_asian <- as.numeric(final1$nh_asian)
library(dplyr)
final2<-subset(final1, select=c(WI,percentfemale,mage,hispan,nh_white,nh_black,nh_asian,SVI,group))%>%
group_by(group)
names(final2)
## [1] "WI" "percentfemale" "mage" "hispan"
## [5] "nh_white" "nh_black" "nh_asian" "SVI"
## [9] "group"
library(data.table)
##
## Attaching package: 'data.table'
## The following objects are masked from 'package:dplyr':
##
## between, first, last
setnames(final2, old = c('percentfemale','mage','hispan','nh_white','nh_black','nh_asian'), new = c('%Female','Median Age','%Hispanics','%Non-hispanic White','%Non-hispanic Black','%Non-hispanic Asian'))
final2
## # A tibble: 5,222 × 9
## # Groups: group [4]
## WI `%Female` `Median Age` `%Hispanics` `%Non-hispanic White`
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 6.83 51.7 21.3 20.3 66.9
## 2 3.33 52.6 39.2 55.4 43.7
## 3 7.67 49.1 37.4 29.7 66.2
## 4 4 50.1 52 15.4 80.6
## 5 8.17 47.2 40.5 42.6 46
## 6 9 52.4 31.7 58.9 30.5
## 7 14.2 50.9 42 29.6 64.4
## 8 7.33 48.7 34.3 94.9 2.9
## 9 4.67 52.0 31.7 90.1 7.1
## 10 7.83 51.3 31.5 23 44
## # ℹ 5,212 more rows
## # ℹ 4 more variables: `%Non-hispanic Black` <dbl>, `%Non-hispanic Asian` <dbl>,
## # SVI <dbl>, group <fct>
final2 <- na.omit(final2)
library(arsenal)
## Warning: package 'arsenal' was built under R version 4.4.2
table_1 <- tableby(group ~ ., data =final2)
summary(table_1, title = "Descriptive Analysis I")
##
##
## Table: Descriptive Analysis I
##
## | | Q1 (N=1323) | Q2 (N=1343) | Q3 (N=1244) | Q4 (N=1297) | Total (N=5207) | p value|
## |:---------------------------|:---------------:|:---------------:|:---------------:|:---------------:|:---------------:|-------:|
## |**WI** | | | | | | < 0.001|
## | Mean (SD) | 4.276 (1.172) | 7.222 (0.693) | 9.504 (0.657) | 12.629 (1.419) | 8.365 (3.241) | |
## | Range | 1.000 - 6.000 | 6.000 - 8.333 | 8.500 - 10.667 | 10.833 - 19.000 | 1.000 - 19.000 | |
## |**%Female** | | | | | | < 0.001|
## | Mean (SD) | 49.655 (5.656) | 50.485 (4.187) | 50.807 (3.932) | 50.390 (4.001) | 50.327 (4.527) | |
## | Range | 1.145 - 66.900 | 4.658 - 63.081 | 17.187 - 71.375 | 26.809 - 68.496 | 1.145 - 71.375 | |
## |**Median Age** | | | | | | < 0.001|
## | Mean (SD) | 38.300 (7.832) | 35.766 (6.341) | 34.846 (5.931) | 35.514 (6.175) | 36.127 (6.749) | |
## | Range | 18.800 - 73.700 | 19.800 - 60.300 | 19.600 - 64.700 | 19.800 - 61.700 | 18.800 - 73.700 | |
## |**%Hispanics** | | | | | | < 0.001|
## | Mean (SD) | 33.477 (28.184) | 38.632 (26.977) | 44.939 (28.221) | 41.626 (27.782) | 39.575 (28.093) | |
## | Range | 0.000 - 100.000 | 0.600 - 100.000 | 0.000 - 100.000 | 0.000 - 99.900 | 0.000 - 100.000 | |
## |**%Non-hispanic White** | | | | | | < 0.001|
## | Mean (SD) | 54.310 (28.583) | 42.348 (26.484) | 34.523 (25.301) | 38.847 (25.518) | 42.646 (27.517) | |
## | Range | 0.000 - 100.000 | 0.000 - 94.400 | 0.000 - 89.800 | 0.000 - 98.400 | 0.000 - 100.000 | |
## |**%Non-hispanic Black** | | | | | | < 0.001|
## | Mean (SD) | 8.646 (13.445) | 12.925 (16.526) | 14.036 (17.738) | 11.210 (13.402) | 11.676 (15.500) | |
## | Range | 0.000 - 95.600 | 0.000 - 96.800 | 0.000 - 92.100 | 0.000 - 81.300 | 0.000 - 96.800 | |
## |**%Non-hispanic Asian** | | | | | | < 0.001|
## | Mean (SD) | 1.502 (3.414) | 3.955 (6.992) | 4.520 (7.157) | 6.031 (8.701) | 3.984 (7.023) | |
## | Range | 0.000 - 42.800 | 0.000 - 61.200 | 0.000 - 54.700 | 0.000 - 75.800 | 0.000 - 75.800 | |
## |**SVI** | | | | | | < 0.001|
## | Mean (SD) | 0.470 (0.254) | 0.504 (0.288) | 0.554 (0.296) | 0.475 (0.308) | 0.500 (0.289) | |
## | Range | 0.000 - 0.996 | 0.003 - 0.999 | 0.002 - 1.000 | 0.000 - 1.000 | 0.000 - 1.000 | |
final3<-subset(combine1, select=c(OBESITY_CrudePrev,LPA_CrudePrev,group))
names(final3)
## [1] "OBESITY_CrudePrev" "LPA_CrudePrev" "group"
final3 <- na.omit(final3)
options(digits=2)
library(arsenal)
table_2 <- tableby(group ~ ., data =final3)
summary(table_2)
##
##
## | | Q1 (N=1332) | Q2 (N=1346) | Q3 (N=1246) | Q4 (N=1298) | Total (N=5222) | p value|
## |:---------------------------|:---------------:|:---------------:|:---------------:|:---------------:|:---------------:|-------:|
## |**OBESITY_CrudePrev** | | | | | | < 0.001|
## | Mean (SD) | 37.376 (4.641) | 36.907 (5.820) | 37.690 (6.372) | 35.365 (6.778) | 36.830 (6.009) | |
## | Range | 19.100 - 52.700 | 22.100 - 52.900 | 21.600 - 55.300 | 18.300 - 53.100 | 18.300 - 55.300 | |
## |**LPA_CrudePrev** | | | | | | < 0.001|
## | Mean (SD) | 27.345 (6.641) | 27.133 (8.374) | 28.904 (8.882) | 26.167 (9.695) | 27.369 (8.509) | |
## | Range | 11.500 - 51.900 | 11.500 - 52.200 | 10.900 - 51.600 | 10.100 - 56.300 | 10.100 - 56.300 | |
Histogram
library(ggplot2)
ggplot(combine1, aes(x = WI, fill = group)) +
geom_histogram(color = "black", bins = 30) +
scale_fill_manual(
values = c("#d9f0d3", "#a6dba0", "#5aae61", "#1b7837"),
name = "Group",
labels = c("Q1", "Q2", "Q3", "Q4")
) +
theme_minimal() +
theme(
legend.position = "bottom",
legend.direction = "horizontal",
legend.title = element_text(size = 12),
legend.text = element_text(size = 10),
legend.spacing.x = unit(1, 'cm'), # <-- Increase space between legend keys
axis.title.x = element_text(size = 14),
axis.title.y = element_text(size = 14),
axis.text = element_text(size = 10),
plot.title = element_text(size = 16, face = "bold")
) +
labs(
x = "Walkability Index",
y = "Number of Census Tracts",
title = "Distribution of Walkability Index by Quartiles"
)

Scatterplots
library(ggplot2)
P1 <- ggplot(combine1a, aes(WI, Obesity)) +
geom_point() +
geom_smooth(method = 'lm') +
labs(x = "Walkability Index", y = "Obesity") +
theme(
axis.title.x = element_text(size = 12.5), # Axis title size
axis.title.y = element_text(size = 12.5),
axis.text.x = element_text(size = 11), # Axis number (tick) size
axis.text.y = element_text(size = 11)
)
P1
## `geom_smooth()` using formula = 'y ~ x'

library(ggplot2)
P2 <- ggplot(combine1a, aes(WI, LPA)) +
geom_point() +
geom_smooth(method = 'lm') +
labs(x = "Walkability Index", y = "Lack of Physical Activity") +
theme(
axis.title.x = element_text(size = 12.5), # Axis title size
axis.title.y = element_text(size = 12.5),
axis.text.x = element_text(size = 11), # Axis number (tick) size
axis.text.y = element_text(size = 11)
)
P2
## `geom_smooth()` using formula = 'y ~ x'

library(gridExtra)
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
plot_list <- list(P1,P2)
do.call("grid.arrange", c(plot_list, ncol = 1))
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'

Spearman’s Correlation Test
cor.test(combine1$OBESITY_CrudePrev, combine1$WI, method = 'spearman',exact=FALSE)
##
## Spearman's rank correlation rho
##
## data: combine1$OBESITY_CrudePrev and combine1$WI
## S = 3e+10, p-value <2e-16
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho
## -0.11
cor.test(combine1$LPA_CrudePrev, combine1$WI, method = 'spearman',exact=FALSE)
##
## Spearman's rank correlation rho
##
## data: combine1$LPA_CrudePrev and combine1$WI
## S = 3e+10, p-value = 6e-06
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho
## -0.063
Regression models
library(broom)
fit<-lm(OBESITY_CrudePrev~WI, data=combine1)
summary(fit)
##
## Call:
## lm(formula = OBESITY_CrudePrev ~ WI, data = combine1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.86 -4.18 -0.44 4.04 18.63
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 38.9129 0.2278 170.8 <2e-16 ***
## WI -0.2492 0.0254 -9.8 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6 on 5220 degrees of freedom
## Multiple R-squared: 0.0181, Adjusted R-squared: 0.0179
## F-statistic: 96.1 on 1 and 5220 DF, p-value: <2e-16
library(broom)
fit<-lm(LPA_CrudePrev~WI, data=combine1)
summary(fit)
##
## Call:
## lm(formula = LPA_CrudePrev ~ WI, data = combine1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.46 -6.54 -1.01 6.00 29.57
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28.7148 0.3250 88.36 < 2e-16 ***
## WI -0.1610 0.0362 -4.44 9.2e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.5 on 5220 degrees of freedom
## Multiple R-squared: 0.00376, Adjusted R-squared: 0.00357
## F-statistic: 19.7 on 1 and 5220 DF, p-value: 9.16e-06
Multivariate Analysis
library(dplyr)
combine2<-combine1%>%
mutate(percentfemale=(fpop*100)/totpop)
names(combine2)
## [1] "TractFIPS" "WI" "OBESITY_CrudePrev"
## [4] "LPA_CrudePrev" "totpop" "fpop"
## [7] "mage" "hispan" "nh_white"
## [10] "nh_black" "nh_asian" "SVI"
## [13] "group" "percentfemale"
combine2$mage <- as.numeric(combine2$mage)
combine2$hispan <- as.numeric(combine2$hispan)
combine2$nh_white <- as.numeric(combine2$nh_white)
combine2$nh_black <- as.numeric(combine2$nh_black)
combine2$nh_asian <- as.numeric(combine2$nh_asian)
library(broom)
fitA<-lm(OBESITY_CrudePrev~WI+mage+percentfemale+SVI+hispan+nh_white+nh_black, data=combine2)
summary(fitA)
##
## Call:
## lm(formula = OBESITY_CrudePrev ~ WI + mage + percentfemale +
## SVI + hispan + nh_white + nh_black, data = combine2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.094 -1.719 0.072 1.804 9.583
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13.18849 0.68368 19.29 <2e-16 ***
## WI -0.23721 0.01235 -19.21 <2e-16 ***
## mage 0.07765 0.00686 11.32 <2e-16 ***
## percentfemale -0.01449 0.00834 -1.74 0.082 .
## SVI 9.83863 0.20824 47.25 <2e-16 ***
## hispan 0.22225 0.00568 39.10 <2e-16 ***
## nh_white 0.15197 0.00575 26.43 <2e-16 ***
## nh_black 0.28808 0.00638 45.17 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.7 on 5199 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.797, Adjusted R-squared: 0.797
## F-statistic: 2.92e+03 on 7 and 5199 DF, p-value: <2e-16
library(broom)
fitB<-lm(LPA_CrudePrev~WI+mage+percentfemale+SVI+hispan+nh_white+nh_black, data=combine2)
summary(fitB)
##
## Call:
## lm(formula = LPA_CrudePrev ~ WI + mage + percentfemale + SVI +
## hispan + nh_white + nh_black, data = combine2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.771 -2.239 -0.037 2.181 15.633
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.22342 0.89626 4.71 2.5e-06 ***
## WI -0.16318 0.01618 -10.08 < 2e-16 ***
## mage 0.13163 0.00899 14.63 < 2e-16 ***
## percentfemale -0.01347 0.01093 -1.23 0.22
## SVI 20.83245 0.27298 76.31 < 2e-16 ***
## hispan 0.14229 0.00745 19.10 < 2e-16 ***
## nh_white 0.06440 0.00754 8.54 < 2e-16 ***
## nh_black 0.14044 0.00836 16.80 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.5 on 5199 degrees of freedom
## (15 observations deleted due to missingness)
## Multiple R-squared: 0.827, Adjusted R-squared: 0.827
## F-statistic: 3.55e+03 on 7 and 5199 DF, p-value: <2e-16
Mediation Analysis:
library(broom)
model1<-lm(LPA_CrudePrev~WI, data=combine1)
summary(model1)
##
## Call:
## lm(formula = LPA_CrudePrev ~ WI, data = combine1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -16.46 -6.54 -1.01 6.00 29.57
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28.7148 0.3250 88.36 < 2e-16 ***
## WI -0.1610 0.0362 -4.44 9.2e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.5 on 5220 degrees of freedom
## Multiple R-squared: 0.00376, Adjusted R-squared: 0.00357
## F-statistic: 19.7 on 1 and 5220 DF, p-value: 9.16e-06
library(broom)
model2<-lm(OBESITY_CrudePrev~WI, data=combine1)
summary(model2)
##
## Call:
## lm(formula = OBESITY_CrudePrev ~ WI, data = combine1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.86 -4.18 -0.44 4.04 18.63
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 38.9129 0.2278 170.8 <2e-16 ***
## WI -0.2492 0.0254 -9.8 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6 on 5220 degrees of freedom
## Multiple R-squared: 0.0181, Adjusted R-squared: 0.0179
## F-statistic: 96.1 on 1 and 5220 DF, p-value: <2e-16
library(broom)
model3<-lm(OBESITY_CrudePrev~LPA_CrudePrev, data=combine1)
summary(model3)
##
## Call:
## lm(formula = OBESITY_CrudePrev ~ LPA_CrudePrev, data = combine1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.971 -1.759 -0.195 1.571 9.135
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 19.45791 0.12274 159 <2e-16 ***
## LPA_CrudePrev 0.63474 0.00428 148 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.6 on 5220 degrees of freedom
## Multiple R-squared: 0.808, Adjusted R-squared: 0.808
## F-statistic: 2.2e+04 on 1 and 5220 DF, p-value: <2e-16
library(psych)
## Warning: package 'psych' was built under R version 4.4.2
##
## Attaching package: 'psych'
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
modelA <- mediate(Obesity ~WI+ (LPA), data =combine1a)

plot(modelA)
MAPPING
library(dplyr)
library(tigris)
## To enable caching of data, set `options(tigris_use_cache = TRUE)`
## in your R script or .Rprofile.
library(leaflet)
tx_censustracts <- tracts(state = "TX")
## Retrieving data for the year 2022
tx_censustracts_cb <- tracts("TX", cb = TRUE)
## Retrieving data for the year 2022
plot(tx_censustracts$geometry)

plot(tx_censustracts_cb$geometry)

dmap <- tx_censustracts_cb %>%
mutate(GEOID = as.numeric(GEOID))
library(dplyr) # Make sure dplyr is loaded
dmap1 <- dmap %>%
left_join(combine1a, by = c("GEOID" = "TractFIPS")) %>%
dplyr::select(GEOID, WI, Obesity, LPA)
names(dmap1)
## [1] "GEOID" "WI" "Obesity" "LPA" "geometry"
library(tmap)
## Breaking News: tmap 3.x is retiring. Please test v4, e.g. with
## remotes::install_github('r-tmap/tmap')
m1 <- tm_shape(dmap1) +
tm_fill(
col = "Obesity",
style = "quantile",
n = 4,
palette = "Reds",
colorNA = "lightgray",
textNA = "Missing",
title = "% Obesity"
) +
tm_borders(col = "gray50", lwd = 0.5) +
tm_scale_bar(
position = c(0.78, 0.04),
breaks = c(0, 100, 200)
) +
tm_compass(
type = "4star",
size = 2.5,
position = c(0.8, 0.12)
) +
tm_layout(
title = "Obesity",
title.position = c(0.65, 0.95),
title.size = 1.2,
title.fontface = "bold",
frame = TRUE,
legend.position = c("left", "bottom"),
legend.title.size = 0.9,
legend.text.size = 0.8
)
m1

library(tmap)
m2 <- tm_shape(dmap1) +
tm_fill(
col = "LPA",
style = "quantile",
n = 4,
palette = "Blues",
colorNA = "lightgray",
textNA = "Missing",
title = "% Lack of Physical Activity"
) +
tm_borders(col = "gray50", lwd = 0.5) +
tm_scale_bar(
position = c(0.78, 0.04), # Moved scale bar more to the right
breaks = c(0, 100, 200)
) +
tm_compass(
type = "4star",
size = 2.5, # Bigger compass
position = c(0.8, 0.12) # Slightly higher than scale bar, more spacing
) +
tm_layout(
title = "Lack of Physical Activity",
title.position = c(0.52, 0.95),
title.size = 1.1,
title.fontface = "bold",
frame = TRUE,
legend.position = c("left", "bottom"),
legend.title.size = 0.9,
legend.text.size = 0.8
)
m2

library(tmap)
m3 <- tm_shape(dmap1) +
tm_fill(
col = "WI",
style = "quantile",
n = 4,
palette = "Greens",
colorNA = "lightgray",
textNA = "Missing",
title = "Walkability Index"
) +
tm_borders(col = "gray50", lwd = 0.5) +
tm_scale_bar(
position = c(0.78, 0.04), # Moved scale bar more to the right
breaks = c(0, 100, 200)
) +
tm_compass(
type = "4star",
size = 2.5, # Bigger compass
position = c(0.8, 0.12) # Slightly higher than scale bar, more spacing
) +
tm_layout(
title = "Walkability",
title.position = c(0.65, 0.95),
title.size = 1.2,
title.fontface = "bold",
frame = TRUE,
legend.position = c("left", "bottom"),
legend.title.size = 0.9,
legend.text.size = 0.8
)
m3
