#Load tidyverse package
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
##
## The following object is masked from 'package:dplyr':
##
## recode
##
## The following object is masked from 'package:purrr':
##
## some
library(knitr)
# read obgyn access data
data <- read.csv("obgyn_access_data.csv")
# Summarize the data
summary(data)
## State.Territory M.D..Obstetrics.and.Gynecology..County.Level.File.
## Length:51 Min. : 56.0
## Class :character 1st Qu.: 214.8
## Mode :character Median : 582.0
## Mean : 799.0
## 3rd Qu.:1025.5
## Max. :4928.0
## NA's :1
## Population.All..County.Level.File. Rate..per.100000.population.
## Min. : 582328 Min. : 6.86
## 1st Qu.: 1854573 1st Qu.:10.07
## Median : 4561284 Median :11.15
## Mean : 6575426 Mean :11.78
## 3rd Qu.: 7625559 3rd Qu.:12.91
## Max. :39368078 Max. :19.20
## NA's :1 NA's :1
## MMR..per.100000.live.births. Uninsured....women.age.15.44.
## Min. :10.10 Min. : 3.000
## 1st Qu.:18.40 1st Qu.: 7.912
## Median :25.40 Median :10.012
## Mean :25.23 Mean :11.015
## 3rd Qu.:30.20 3rd Qu.:13.506
## Max. :43.50 Max. :24.550
## NA's :6 NA's :7
## Median.Household.Income High.School.Completion.. X..of.C.sections
## Min. :45928 Min. :83.30 Min. :22.77
## 1st Qu.:57343 1st Qu.:87.05 1st Qu.:28.07
## Median :62028 Median :89.85 Median :30.85
## Mean :64738 Mean :89.20 Mean :30.72
## 3rd Qu.:74823 3rd Qu.:91.33 3rd Qu.:33.71
## Max. :86644 Max. :93.60 Max. :38.38
## NA's :7 NA's :7 NA's :7
## Late.No.Prenatal.Care.. black_mmr
## Min. : 3.450 Min. : 3.60
## 1st Qu.: 4.844 1st Qu.: 30.85
## Median : 5.975 Median : 54.90
## Mean : 6.253 Mean : 54.96
## 3rd Qu.: 7.088 3rd Qu.: 74.15
## Max. :11.750 Max. :122.90
## NA's :7 NA's :1
Since there are some NA’s present and the column labels are lengthy, the data will be recoded
# recode the data
data.clean <- data %>%
rename(state = State.Territory) %>%
rename(total_obgyn = M.D..Obstetrics.and.Gynecology..County.Level.File.) %>%
rename(population = Population.All..County.Level.File.) %>%
rename(obgyn_per_100K = Rate..per.100000.population.) %>%
rename(mmr = MMR..per.100000.live.births.) %>%
rename(unins = Uninsured....women.age.15.44.) %>%
rename(income = Median.Household.Income) %>%
rename(education = High.School.Completion..) %>%
rename(c.sect = X..of.C.sections) %>%
rename(prenatal = Late.No.Prenatal.Care..) %>%
drop_na()
# Repeat data summary
summary(data.clean)
## state total_obgyn population obgyn_per_100K
## Length:44 Min. : 77.0 Min. : 731158 Min. : 6.86
## Class :character 1st Qu.: 277.8 1st Qu.: 3014588 1st Qu.:10.03
## Mode :character Median : 673.5 Median : 5437691 Median :11.15
## Mean : 890.4 Mean : 7336484 Mean :11.60
## 3rd Qu.:1151.8 3rd Qu.: 8663515 3rd Qu.:12.53
## Max. :4928.0 Max. :39368078 Max. :19.20
## mmr unins income education
## Min. :10.10 Min. : 3.000 Min. :45928 Min. :83.30
## 1st Qu.:18.38 1st Qu.: 7.912 1st Qu.:57343 1st Qu.:87.05
## Median :25.55 Median :10.012 Median :62028 Median :89.85
## Mean :25.27 Mean :11.015 Mean :64738 Mean :89.20
## 3rd Qu.:30.23 3rd Qu.:13.506 3rd Qu.:74823 3rd Qu.:91.33
## Max. :43.50 Max. :24.550 Max. :86644 Max. :93.60
## c.sect prenatal black_mmr
## Min. :22.77 Min. : 3.450 Min. : 3.60
## 1st Qu.:28.07 1st Qu.: 4.844 1st Qu.: 33.77
## Median :30.85 Median : 5.975 Median : 56.50
## Mean :30.72 Mean : 6.253 Mean : 58.08
## 3rd Qu.:33.71 3rd Qu.: 7.088 3rd Qu.: 74.92
## Max. :38.38 Max. :11.750 Max. :122.90
The data is much more polished and managable now
# Plot a histogram of obgyn access level
ggplot(data = data.clean, aes(x = obgyn_per_100K)) +
geom_histogram(bins = 15)
Interpretation: Number of OBGYNs per 100K population is not normally
distributed
# Plot histogram of maternal mortality rate
data.clean %>%
ggplot(aes(x = mmr)) +
geom_histogram(bins = 15)
Interpretation: Maternal mortality rate is not normally distributed
across states
# Create a scatter plot of OBGYNs per 100K population and maternal mortality rate with a linear fit line
data.clean %>%
ggplot(aes(x = obgyn_per_100K, y = mmr)) +
geom_point(aes(size = "state"), color = "blue") +
geom_smooth(aes(linetype = "Linear fit line"), method = "lm", se = FALSE, color = "gray60") +
theme_minimal() +
labs(x = "Number of OBGYNs per 100K Population", y = "Maternal Mortality Rate") +
ggtitle("Maternal mortality rate vs. OBGYNs per 100K")
## Warning: Using size for a discrete variable is not advised.
## `geom_smooth()` using formula = 'y ~ x'
Interpretation: The is a weakly negative correlation between # of OBGYNs
per 100K population and MMR.
# Create a linear model comparing MMR to OBGYN access
prd.by.ob <- lm(formula = mmr ~ obgyn_per_100K, data = data.clean, na.action = na.exclude)
# Get linear model summary
summary(prd.by.ob)
##
## Call:
## lm(formula = mmr ~ obgyn_per_100K, data = data.clean, na.action = na.exclude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -14.7271 -6.1656 0.2045 3.9359 16.9358
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 34.6681 5.9911 5.787 8.05e-07 ***
## obgyn_per_100K -0.8098 0.5044 -1.605 0.116
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.497 on 42 degrees of freedom
## Multiple R-squared: 0.05781, Adjusted R-squared: 0.03538
## F-statistic: 2.577 on 1 and 42 DF, p-value: 0.1159
The y intercept is 34.6681. This estimates that a state with 0 obgyns would have 34.6681 pregnancy related deaths per 100K births.
The slope is -0.8098. This means if the number of obgyns per 100K population in a state goes up by 1, there would be 0.8098 fewer pregnancy related deaths per 100K births.
# Get the confidence intervals
confint(object = prd.by.ob)
## 2.5 % 97.5 %
## (Intercept) 22.577608 46.7585049
## obgyn_per_100K -1.827799 0.2082004
Interpretation: A simple linear regression analysis found that a state’s number of OBGYNs is not a statistically significant predictor for number of pregnancy related deaths (b = -0.8098, p > 0.05). The slope value is likely between -1.828 and 0.208. The model was not statistically better from the baseline (F(1, 42) = 2.577, p > 0.05). These results do not suggest a relationship between number obgyns and number of pregnancy related deaths in a state.
Checking assumptions
1: Independent Observations * Number of pregancy related deaths are not reliant on any other pregnancy related death. * This assumption is met
2: Continuous outcome * Since NAs were eliminated, # of OBGYNs and # of pregnancy related deaths can be any value greater than or equal to 0 * This assumption is met
3: Linearity
# Plot a scatter plot of obgyn access vs mmr with a linear fit line and loess curve
data.clean %>%
ggplot(aes(x = obgyn_per_100K, y = mmr)) +
geom_point() +
geom_smooth(aes(color = "Linear fit line"), method = "lm", se = FALSE) +
geom_smooth(aes(color = "Loess curve"), se = FALSE) +
theme_minimal() +
labs(y = "Maternal Mortality Rate", x = "OBGYNs per 100K population") +
scale_color_manual(values = c("gray60", "deeppink"), name = "") +
ggtitle("Maternal mortality rate vs. OBGYNs per 100K population")
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Interpretation: The Loess curve does not match the linear fit line at all. The assumption is failed.
4: Homeoscadacity
# Conduct a bp test for homeoscadacity
const.var.test <- lmtest::bptest(formula = prd.by.ob)
const.var.test
##
## studentized Breusch-Pagan test
##
## data: prd.by.ob
## BP = 0.4322, df = 1, p-value = 0.5109
The Breush-Pagan test has a high p-value. The null hypothesis is accepted, and the assumption of homoescedasticity is met.
Assumption 5: Independence of residuals
# Perforn a Durbin-Watson test for independence of residuals
lmtest::dwtest(formula = prd.by.ob)
##
## Durbin-Watson test
##
## data: prd.by.ob
## DW = 1.8048, p-value = 0.2533
## alternative hypothesis: true autocorrelation is greater than 0
The D-W statistic is near 2 and the p-value is > 0.05, so the null hypothesis is accepted. This assumption is met.
Assumption 6: Normality of residuals
# Plot a histogram of the residuals
data.frame(prd.by.ob$residuals) %>%
ggplot(aes(x = prd.by.ob.residuals)) +
labs(x = "Residuals") +
geom_histogram(bins = 15) +
ggtitle("Residuals of MMR vs. OBGYNs per 100K")
The residuals are not normally distributed. This assumption is not
met.
Conclusion: A linear model is not an appropriate way to describe this data
Since a linear fit model is not an appropriate choice for analysis, covariance will be tried instead.
Check the strength of the Pearson’s Product Moment Correlation Coefficient
# get the Pearson's product moment correlation coefficient
data.clean %>% summarize(cov.ob.prd = cor(x = obgyn_per_100K,
y = mmr,
use = "complete"))
## cov.ob.prd
## 1 -0.2404432
Interpretation: the strength of the r is weakly negative
Check assumptions (see above for tests)
1: Observations are independent: met 2: Both variables are continuous: met 3: Both variables are normally distributed: not met 4: Linearity: not met 5: Homoscedasticity: met
Because one of the the main issue is the nonlinaerity of the relationship, I will first try transforming the data
# apply log and arcsin transformations to the data
data.clean <- data.clean %>%
mutate(log.ob = log(x = (obgyn_per_100K/100)/(1-obgyn_per_100K/100))) %>%
mutate(log.prd = log(x = (mmr/100)/(1-mmr/100))) %>%
mutate(arcsin.ob = asin(x = sqrt(obgyn_per_100K/100))) %>%
mutate(arcsin.prd = asin(x = sqrt(mmr/100)))
Recheck normality
# Check normality of log transformed obgyn access using a box plot
data.clean %>%
ggplot(aes(x = log.ob)) +
geom_histogram(binwidth = 0.05)
Interpretation: a log transformation does not make the obgyn access data
more normally distributed
# Check normality of log transformed mmr with a histogram
data.clean %>%
ggplot(aes(x = log.prd)) +
geom_histogram(bins = 15)
Interpretation: Log transformation was unsuccessful in normalizing mmr
data
# Check normality of arcsin transformed obgyn access data with a histogram
data.clean %>%
ggplot(aes(x = arcsin.ob)) +
geom_histogram(bins = 15)
Interpretation: Arsin transformation was unsuccessful in normalizing
obgyn access data
# Check the normality of arcsin transformed mmr data with a histogram
data.clean %>%
ggplot(aes(x = arcsin.prd)) +
geom_histogram(bins = 15)
Interpretation: Arsin transformation was unsuccessful in normalizing mmr data
Next, folded power transformation will be attempted
# Get lambda tranformation value for obgyn access
p.ob <- rcompanion::transformTukey(x = data.clean$obgyn_per_100K,
plotit = FALSE,
quiet = TRUE,
returnLambda = TRUE)
p.ob
## lambda
## -0.25
# Get lambda tranformation value for mmr
p.prd <- rcompanion::transformTukey(x = data.clean$mmr,
plotit = FALSE,
quiet = TRUE,
returnLambda = TRUE)
p.prd
## lambda
## 0.275
# Apply folded power transformation to data
data.clean <- data.clean %>%
mutate(folded.ob = (obgyn_per_100K/100)^(1/p.ob)-(1-obgyn_per_100K/100)^(1/p.ob)) %>%
mutate(folded.prd = (mmr/100)^(1/p.prd)-(1-mmr/100)^(1/p.prd))
# Check normality of folded power transformed obgyn access data for normality
data.clean %>%
ggplot(aes(x = folded.ob)) +
geom_histogram(bins = 15)
Interpretation: Folded power transformation was unsuccessful at
normalizing obgyn access data
# Check normality of folded power transformed mmr data using a histogram
data.clean %>%
ggplot(aes(x = folded.prd)) +
geom_histogram(bins = 15)
Interpretation: Folded power transformation was unsuccessful at normalizing mmr data
Conclusion: Transforming the data was not successful. It may be better to transform the data into categories.
Recoding data as categorical
OBGYN access data will be recoded as categorical, high OBGYN access, mid OBGYN access, and low OBGYN access. I will do this by getting the range of the data, and dividing by 3. This will provide the low, mid, and high cutoffs.
# Define the cutoff value for
cutoff <- (max(data.clean$obgyn_per_100K) - min(data.clean$obgyn_per_100K))/3;cutoff
## [1] 4.113333
# Define the minimum values for each obgyn access category
low.ob <- min(data.clean$obgyn_per_100K);low.ob
## [1] 6.86
mid.ob <- low.ob + cutoff;mid.ob
## [1] 10.97333
high.ob <- mid.ob + cutoff;high.ob
## [1] 15.08667
# Add a column to the dataset for obgyn access category
data.clean <- data.clean %>%
mutate(ob.category = 0)
# Population the newly created ob.category variable based on OBGYNs per 100K
for(i in 1:nrow(data.clean)){
if(data.clean$obgyn_per_100K[i] < mid.ob){
data.clean$ob.category[i] <- "low"}
else if (data.clean$obgyn_per_100K[i] >= mid.ob & data.clean$obgyn_per_100K[i] < high.ob){
data.clean$ob.category[i] <- "mid"
} else{
data.clean$ob.category[i] <- "high"
}
}
The appropriate parametric test for this data is ANOVA. First, assumptions will be tested.
1: Independence of observations: met (see above) 2: Continuous variables and independent groups: met 3: Normal distribution in each group
# Perform Shapiro Wilk test for normality among groups
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = mmr)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <chr> <dbl>
## 1 high 0.0367
## 2 low 0.563
## 3 mid 0.569
For high obgyn access states, mmr is not normally distributed (Shapiro-Wilk test, p < 0.05). This assumption is failed.
4: Equal variances for each group
# Perform a Levene test for homogeny of variance
leveneTest(mmr ~ ob.category, data.clean)
## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.2444 0.7843
## 41
A Levene test confirms variance is equal between groups (p >0.05).
Since the data fails the normality assumption, we will use the Kruskal-Wallis rank sum test.
# Run KW test comparing mmr to obgyn access category
kw.ob.by.prd <- kruskal.test(formula = mmr ~ ob.category, data = data.clean);kw.ob.by.prd
##
## Kruskal-Wallis rank sum test
##
## data: mmr by ob.category
## Kruskal-Wallis chi-squared = 2.8224, df = 2, p-value = 0.2438
Based on a Kruskal-Wallis test, we accept the null hypothesis that there is no difference in maternal mortality between states with low, mid, and high obgyn availability (X2 = 2.8224, p > 0.05).
# order the obgyn access categories as factors
data.clean$ob.category <- factor(data.clean$ob.category, levels = c("low", "mid", "high"))
# Create a box plot of mmr in states with low, mid, and high OBGYN access
data.clean %>%
ggplot(aes(x = ob.category, y = mmr)) +
geom_boxplot(aes(fill = ob.category)) +
theme_minimal() +
scale_fill_brewer(palette = "Spectral", guide = "none") +
scale_color_brewer(palette = "Spectral", guide = "none") +
labs(x = "Level of obgyn access", y = "Maternal mortality rate") +
ggtitle("Maternal mortality rate v. OBGYN access level")
Next, social determininants of health will be compared to OBGYN access level. Assumptions for ANOVA will be considered first, then an appropriate test will be determined.
How does percent of women age 15-44 in a state effect obgyn access?
1: Independece of observations: met
2: Continuous variable and independent groups: Met; percent uninsured can range from 0-100 continuously. OBGYN access levels are independent of other states.
3: Assumption of normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = unins)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.758
## 2 mid 0.128
## 3 high 0.668
All categories have p > 0.05 for the Shapiro-Wilk test. The null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
# Perform a Levene test for homogeny of variance
leveneTest(unins ~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.3371 0.2738
## 41
With p > 0.05, the Levene Test confirms that the variance between groups is homogenous. The assumption is met.
Based on this data, I will proceed with a Kruskall-Wallis test again
# Create a boxplot of percent of women ages 15-44 uninsured grouped by obgyn access level
data.clean %>%
ggplot(aes(x = ob.category, y = unins)) +
geom_boxplot(aes(fill = ob.category)) +
scale_fill_brewer(palette = "Spectral", guide = "none") +
scale_color_brewer(palette = "Spectral", guide = "none") +
theme_minimal() +
labs(x = "Level of obgyn access", y = "% Women 15-44 uninsured") +
ggtitle("% Women age 15-44 Uninsured vs. OBGYN access level")
#
kw.ob.by.unins <- kruskal.test(formula = unins ~ ob.category, data = data.clean);kw.ob.by.unins
##
## Kruskal-Wallis rank sum test
##
## data: unins by ob.category
## Kruskal-Wallis chi-squared = 7.3577, df = 2, p-value = 0.02525
Based on a Kruskal-Wallis test, the null hypothesis is rejected in favor of the alternate hypothesis that a relationship exists between percent of the population that is uninsured and access to obgyn care (X2 = 7.3577, p < 0.05).
# Run a dunn post hoc test on the relationship between the % population of uninsured women and obgyn access level
dunn.test::dunn.test(x = data.clean$unins, g = data.clean$ob.category, method = "bonferroni")
## Kruskal-Wallis rank sum test
##
## data: x and group
## Kruskal-Wallis chi-squared = 7.3577, df = 2, p-value = 0.03
##
##
## Comparison of x by group
## (Bonferroni)
## Col Mean-|
## Row Mean | high low
## ---------+----------------------
## low | -2.709367
## | 0.0101*
## |
## mid | -2.077258 0.969535
## | 0.0567 0.4984
##
## alpha = 0.05
## Reject Ho if p <= alpha/2
Based on the Dunn test, there is a significant difference in percent of the population uninsured in states with high vs. low obgyn access (p < 0.05).
Effect size
etaSquared.unins.ob <- (7.3577-3+1)/(44+3);etaSquared.unins.ob
## [1] 0.1139936
The effect size is medium
A Kruskal-Wallis test found a statistically significant difference in percent of the population uninsured across states with different levels of obgyn access (X2 = 7.3577, p < 0.05). Based on a Dunn’s post hoc test, states with low obgyn access had a significantly higher percent of the population uninsured than states with high obgyn access (p < 0.05). There was a medium effect size for the relationship between obgyn access level and ranked values of population uninsured (n2 = 0.114).
Does OBGYN access differ by median income
Check assumptions
1: Independent observations: met; one household’s income does not affect another’s 2: Continuous variable and independent groups: met; median household income can be any positive value. OBGYN access level is not dependent on access level of other states.
3: Normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = income)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.341
## 2 mid 0.0527
## 3 high 0.258
All groups have p > 0.05 for a shapiro-wilk test, so the null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
leveneTest(income ~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.6546 0.525
## 41
The Levene’s test p-value is > 0.05, so the null hypothesis that variance is homogenous is met.
ANOVA is the appropriate analysis for this data.
H0: There is no relationship between median household income and obgyn access level in a state. HA: A relationship exists between median household income and obgyn access level in a state.
# Perform ANOVA
aov.inc.ob <- aov(income ~ ob.category, data = data.clean)
summary(aov.inc.ob)
## Df Sum Sq Mean Sq F value Pr(>F)
## ob.category 2 1.282e+09 641163427 6.583 0.00332 **
## Residuals 41 3.993e+09 97393894
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Based on an ANOVA test, the null hypothesis is rejected in favor of the alternate hypothesis, that a relationship exists between median household income and obgyn access level in a state (F = 6.583, p < 0.01).
# Make a box plot showing distribution of income grouped by obgyn access
data.clean %>%
ggplot(aes(x = ob.category, y = income)) +
geom_boxplot(aes(fill = ob.category)) +
scale_fill_brewer(palette = "Spectral", guide = "none") +
scale_color_brewer(palette = "Spectral", guide = "none") +
theme_minimal() +
labs(x = "Level of obgyn access", y = "Median Income") +
ggtitle("Median income vs. level of OBGYN access")
As median household income increases, so does obgyn access.
A Bonferroni pairwise test will be run as a post hoc test
# Run Bonferroni pairwise post hoc test.
pairwise.t.test(x = data.clean$income,
g = data.clean$ob.category,
p.adj = "bonf")
##
## Pairwise comparisons using t tests with pooled SD
##
## data: data.clean$income and data.clean$ob.category
##
## low mid
## mid 0.025 -
## high 0.011 0.589
##
## P value adjustment method: bonferroni
Based on the Bonferroni pairwise comparison, states with low obgyn access have lower median household income than states with high obgyn access (p < 0.05) and states with mid obgyn access (p < 0.05).
Effect size
# Calculate omega squared effect size
omega.inc.ob <- (6.583-1)/(6.583 + ((44-3+1)/3-1));omega.inc.ob
## [1] 0.2850942
Based on the omega squared value, the effect size is large (w2 = 0.285).
A one-way ANOVA found a statistically significant difference in percent of the population uninsured across states with different levels of obgyn access (F = 6.583, p < 0.01). Based on a Bonferroni post hoc test, states with low obgyn access have significantly lower median income than states with high obgyn access (p < 0.05) and states with mid obgyn access (p < 0.05). There was a large effect size for the relationship between obgyn access level and ranked values of median income (w2 = .285).
Is there are realationship between percent of deliveries done by c-section and obgyn access
1: Independent observations: met; one person’s delivery does not impact another’s 2: Continuous variable and independent groups: met; percent of deliveries by c-section can be any value 0-100. OBGYN access level is not dependent on access level of other states.
3: Normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = c.sect)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.598
## 2 mid 0.267
## 3 high 0.182
All groups have p > 0.05 for a shapiro-wilk test, so the null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
# Perform a Levene test for homogeneity of variance
leveneTest(c.sect ~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 2.3739 0.1058
## 41
The Levene’s test p-value is > 0.05, so the null hypothesis that variance is homogenous is met.
ANOVA is the appropriate analysis for this data.
H0: There is no relationship between percent of births caried out by c-section and obgyn access level in a state. HA: A relationship exists between percent of births caried out by c-section and obgyn access level in a state.
# Perform ANOVA
aov.csect.ob <- aov(c.sect ~ ob.category, data = data.clean)
summary(aov.csect.ob)
## Df Sum Sq Mean Sq F value Pr(>F)
## ob.category 2 37.3 18.65 1.345 0.272
## Residuals 41 568.4 13.86
Based on a one way ANOVA, the null hypothesis is accepted, that no relationship between percent of births caried out by c-section and obgyn access level in a state.
Is there a relationship between percent of population who has completed high school and OBGYN access?
1: Independent observations: met; one person’s education status is not dependent on another’s 2: Continuous variable and independent groups: met; percent of the population that completed high school can be any value 0-100. OBGYN access level is not dependent on access level of other states.
3: Normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = education)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.0793
## 2 mid 0.0683
## 3 high 0.504
All groups have p > 0.05 for a shapiro-wilk test, so the null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
# Perform a Levene test for homogeneity of variance
leveneTest(education ~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.6183 0.5438
## 41
The Levene’s test p-value is > 0.05, so the null hypothesis that variance is homogenous is met.
ANOVA is the appropriate analysis for this data.
H0: There is no relationship between percent of the population with a high school education and obgyn access level in a state. HA: A relationship exists between percent of the population with a high school education and obgyn access level in a state.
# Perform ANOVA
aov.edu.ob <- aov(education ~ ob.category, data = data.clean)
summary(aov.edu.ob)
## Df Sum Sq Mean Sq F value Pr(>F)
## ob.category 2 0.84 0.421 0.055 0.946
## Residuals 41 311.96 7.609
Based on a one way ANOVA, the null hypothesis that no relationship exists between percent of the population with a high school education and obgyn access level in a state is accepted (F = 0.055, p > 0.05).
Is there a relationship between OBGYN access and percent of pregnancies with late/no prenatal care
1: Independent observations: met; one person’s access to or choice to obtain prenatal care is not dependent on another person 2: Continuous variable and independent groups: met; percent of pregnancies with late/no prenatal carel can be any value 0-100. OBGYN access level is not dependent on access level of other states.
3: Normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = prenatal)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.0914
## 2 mid 0.802
## 3 high 0.611
All groups have p > 0.05 for a shapiro-wilk test, so the null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
# Perform a Levene test for homogeneity of variance
leveneTest(prenatal~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.6028 0.5521
## 41
The Levene’s test p-value is > 0.05, so the null hypothesis that variance is homogenous is met.
ANOVA is the appropriate analysis for this data.
H0: There is no relationship between percent of pregnancies with late/no prenatal care and obgyn access level in a state. HA: A relationship exists between percent of pregnancies with late/no prenatal care education and obgyn access level in a state.
# Perform ANOVA
aov.prenat.ob <- aov(prenatal ~ ob.category, data = data.clean)
summary(aov.prenat.ob)
## Df Sum Sq Mean Sq F value Pr(>F)
## ob.category 2 3.34 1.668 0.399 0.673
## Residuals 41 171.30 4.178
Based on one-way analysis of variance, the null hypothesis is accepted, that there is no relationship between percent of pregnancies with late/no prenatal care and obgyn access level in a state (F = 0.399, p > 0.05).
Does OBGYN access impact Black maternal mortality rates?
Check assumptions for ANOVA
1: Independent observations: met; no case of mortality is dependent on another 2: Continuous variable and independent groups: met; MMR can be any positive value. OBGYN access level is not dependent on access level of other states.
3: Normality between groups
# Perform a Shapiro-Wilk test
data.clean %>%
group_by(ob.category) %>%
summarize(shapiro.pval = shapiro.test(x = black_mmr)$p.value)
## # A tibble: 3 × 2
## ob.category shapiro.pval
## <fct> <dbl>
## 1 low 0.906
## 2 mid 0.0918
## 3 high 0.890
All groups have p > 0.05 for a shapiro-wilk test, so the null hypothesis that data is normally distributed between groups is met.
4: Homogeny of variance
# Perform a Levene test for homogeneity of variance
leveneTest(black_mmr ~ ob.category, data.clean)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.091 0.9132
## 41
The Levene’s test p-value is > 0.05, so the null hypothesis that variance is homogenous is met.
ANOVA is the appropriate analysis for this data.
H0: There is no relationship between Black maternal mortality and obgyn access level in a state. HA: A relationship exists between Black maternal mortality and obgyn access level in a state.
# Perform ANOVA
aov.bmmr.ob <- aov(black_mmr ~ ob.category, data = data.clean)
summary(aov.bmmr.ob)
## Df Sum Sq Mean Sq F value Pr(>F)
## ob.category 2 552 275.8 0.322 0.727
## Residuals 41 35166 857.7
Based on one-way analysis of variance, the null hypothesis is accepted, that there is no relationship between Black maternal mortality and obgyn access level in a state (F = 0.322, p > 0.05).