I was working on finding a meaningful metric to explain alcohol access in neighborhoods. One of the assumptions I made was that different types of alcohol outlets have different implications in social characteristics and economic strength of a neighborhood. I am curious to find out whether this assumption is supported by statistics.
Prepare the Dataset
For this task, I am associating each liquor license with a few selected social economic metrics at the census block group level, i.e. assigning each liquor license an value of the indicator based on which census block group it is located; then the metrics are aggregated by census block groups, and finally, I use t-test and ANOVA to analyze the differences in these metrics between each liquor license type.
I have picked median house income and the percentage of white residents as the two indicators for the social fabric of the neighborhoods.
source("src/preparation.R") source("src/alcohol.R") bizz.cols <- c("CT_ID_10", "ALC_BIZ_TYPE", "ALC_TYPE") acs.cols <- c("CT_ID_10", "MedHouseIncome", "White") dat <- merge(bizz.alc[, bizz.cols], acs.ct[, acs.cols], all.x = TRUE)
When categorized by type of alcohol allowed, the licenses can be divided into two group: malt & wine, and all alcoholic beverages. I used a Student’s t-test to see whether census indicators for these two groups are significantly different from each other.
First, the median house income:
library(broom) library(knitr) dat.t <- t.test(MedHouseIncome ~ ALC_TYPE, dat) %>% tidy() dat.t <- dat.t[, c(1:6)] colnames(dat.t) <- c( "mean difference", "mean in group (All Alcohol)", "mean in group (Malte & Wine)", "t statistic", "p value", "df" ) row.names(dat.t) <- c("value") round(dat.t, 2) %>% t() %>% kable()
|mean in group (All Alcohol)||76030.47|
|mean in group (Malte & Wine)||65265.46|
We get a t-test of 5.76, and the p value is as low as 1.1e-8 indicating a very strong difference. On average, a census block group containing full liquor licenses have a $10,765 higher median house income than those with malt & wine only liquor licenses. Note that each group has different number of members (liquor licenses)
This makes sense because full liquor licenses are much harder to get and have higher prices on the market. They naturally would be more distributed among rich neighborhoods.
Then let’s look at race:
dat.t <- t.test(White ~ ALC_TYPE, dat) %>% tidy() dat.t <- dat.t[, c(1:6)] colnames(dat.t) <- c( "mean difference", "mean in group (All Alcohol)", "mean in group (Malte & Wine)", "t statistic", "p value", "df" ) row.names(dat.t) <- c("value") round(dat.t, 2) %>% t() %>% kable()
|mean in group (All Alcohol)||0.66|
|mean in group (Malte & Wine)||0.64|
The p value is 0.22, not so significant. Statistically speaking, race is not related to the distribution of different alcohol types.
We have 6 different business types designated in the categorization of liquor licenses. An ANOVA test can be used to test whether the implications of these business types are statistically significant.
dat.anova <- aov(MedHouseIncome ~ ALC_BIZ_TYPE, dat) broom::tidy(dat.anova) %>% knitr::kable()
With F value’s p value as low as 0.00201, it is statistically significant that at least one business type is different than others.
And based on Tukey’s “Honest Significant Difference” method, three group combinations were identified with significant differences.
broom::tidy(TukeyHSD(dat.anova))[2:6] %>% knitr::kable()
|Farmer Distillery-Common Vectualler||-8965.987||-69573.924||51641.949||0.9982903|
|General On-Premises-Common Vectualler||-14653.944||-32721.074||3413.186||0.1886058|
|General On-Premises-Farmer Distillery||-5687.957||-68810.105||57434.192||0.9998475|
They are Hotel-Club, Hotel-CommonVectualler, and Hotel-GeneralOnPremises.
broom::tidy(TukeyHSD(dat.anova))[c(4, 8, 13), 2:6] %>% knitr::kable()
They are all related to hotels, which indicates that liquor licenses issues to hotels might be the outlier here, i.e., hotels offering on-premise consumption alcohol beverages tend to be located in richer neighborhoods.
The differences between hotel and the remaining two categories, Farmer Distillery and Tavern, are not statistically significant. However, this might just be because of the very small number of licenses in these two categories (with just 2 and 6 each)–the degree of freedom is really low, making it impossible to get a significant result.
Explained by Plots
Using bar chart and error bars, the mean differences of median house income between different liquor licenses can be plotted as following:
library(reshape2) dat.melted <- melt(dat[, c("ALC_BIZ_TYPE", "MedHouseIncome")], id.vars = "ALC_BIZ_TYPE") dat.means <- aggregate(value ~ ALC_BIZ_TYPE + variable, dat.melted, FUN = mean) names(dat.means)<-"mean" dat.ses <- aggregate(value ~ ALC_BIZ_TYPE + variable, dat.melted, function(x) sd(x, na.rm = TRUE) / sqrt(length(!is.na(x))) ) names(dat.ses)<-"se" dat.means <- merge(dat.means, dat.ses, by = c("ALC_BIZ_TYPE", "variable")) dat.means <- transform(dat.means, lower = mean - se, upper = mean + se) ggplot(data = dat.means, aes(x = ALC_BIZ_TYPE, y = mean)) + geom_bar(stat = "identity", position = "dodge", fill = "tomato2") + geom_errorbar( aes(ymax = upper, ymin = lower), position = position_dodge(.9)) + theme(axis.text.x = element_text(angle = 70, hjust = 1)) + ylab("Mean") + xlab("Business Type")
The high income level for liquor licenses in hotels is very obvious.