IPSOS Analysis

Install Packages and Read in Data

#install.packages("foreign")
library(foreign)
library(likert)
library(ggplot2)
library(usmap)
library(tibble)
library(plyr)
library(dplyr)


data <- read.spss("KP_OMNI_2405_BMW_Client_File_02122024.sav", to.data.frame=TRUE)
str(data)

Participants.

part.table=data %>%
  group_by(ppstaten) %>%
  summarize(count=n())
part.table=as.data.frame(part.table)

part.table=part.table %>% rename(state = ppstaten)


plot_usmap(data=part.table, regions = "states", values="count") + 
  labs(title = "U.S. States",
       subtitle = "Proportion of gun owners by state.") + 
  theme(panel.background=element_blank())+
   scale_fill_continuous(low = "antiquewhite1", high ="antiquewhite4", 
                          name = "Number of Participants",
                          limits = c(0,110),
                         na.value="white") +
  theme(legend.position = "right")

BMW1. Do you currently own or have access to a firearm? (All)

summary(data$BMW1)

##           Skipped               Yes                No Prefer not to say 
##                 5               363               643                 6

Q1.comp=subset(data, BMW1 == "Yes" | BMW1 == "No")
summary(Q1.comp$BMW1)

##           Skipped               Yes                No Prefer not to say 
##                 0               363               643                 0

#Gender Subset
chisq.test(Q1.comp$BMW1, Q1.comp$ppgender)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  Q1.comp$BMW1 and Q1.comp$ppgender
## X-squared = 19.556, df = 1, p-value = 9.771e-06

Q1.table <- Q1.comp %>%
  group_by(BMW1, ppgender) %>%
    summarize(count = n())

ggplot(Q1.table, aes(x=BMW1, y=count, fill=ppgender)) + geom_bar(position = "fill", stat="identity")

#party
Q1.party.table=Q1.comp %>%
  group_by(BMW1, xparty4) %>%
  summarize(count=n())
Q1.party.table

## # A tibble: 9 × 3
## # Groups:   BMW1 [2]
##   BMW1  xparty4        count
##   <fct> <fct>          <int>
## 1 Yes   Republican       148
## 2 Yes   Democrat          59
## 3 Yes   Independent      125
## 4 Yes   Something else    31
## 5 No    Republican       119
## 6 No    Democrat         247
## 7 No    Independent      214
## 8 No    Something else    60
## 9 No    Missing            3

ggplot(Q1.party.table, aes(x=BMW1, y=count, fill=xparty4)) + geom_bar(position = "fill", stat="identity")

#urbanicity
Q1.urb.table=Q1.comp %>%
  group_by(BMW1, xurbanicity) %>%
  summarize(count=n())
Q1.urb.table

## # A tibble: 6 × 3
## # Groups:   BMW1 [2]
##   BMW1  xurbanicity count
##   <fct> <fct>       <int>
## 1 Yes   Urban          95
## 2 Yes   Rural         100
## 3 Yes   Suburban      168
## 4 No    Urban         235
## 5 No    Rural          65
## 6 No    Suburban      343

ggplot(Q1.urb.table, aes(x=BMW1, y=count, fill=xurbanicity)) + geom_bar(position = "fill", stat="identity")

#age
ggplot(Q1.comp, aes(x=BMW1, y=ppage)) + geom_boxplot()

#education
Q1.ed.table=Q1.comp %>%
  group_by(BMW1, ppeducat) %>%
  summarize(count=n())
Q1.ed.table

## # A tibble: 8 × 3
## # Groups:   BMW1 [2]
##   BMW1  ppeducat           count
##   <fct> <fct>              <int>
## 1 Yes   Less than HS          16
## 2 Yes   HS                   104
## 3 Yes   Some college         110
## 4 Yes   Bachelor or higher   133
## 5 No    Less than HS          56
## 6 No    HS                   164
## 7 No    Some college         153
## 8 No    Bachelor or higher   270

ggplot(Q1.ed.table, aes(x=BMW1, y=count, fill=ppeducat)) + geom_bar(position = "fill", stat="identity")

#ethnicity
Q1.eth.table=Q1.comp %>%
  group_by(BMW1, ppethm) %>%
  summarize(count=n())
Q1.eth.table

## # A tibble: 10 × 3
## # Groups:   BMW1 [2]
##    BMW1  ppethm                                  count
##    <fct> <fct>                                   <int>
##  1 Yes   White, Non-Hispanic                       287
##  2 Yes   Black or African American, Non-Hispanic    24
##  3 Yes   Other, Non-Hispanic                         9
##  4 Yes   Hispanic                                   28
##  5 Yes   2+ races, Non-Hispanic                     15
##  6 No    White, Non-Hispanic                       409
##  7 No    Black or African American, Non-Hispanic    75
##  8 No    Other, Non-Hispanic                        39
##  9 No    Hispanic                                   95
## 10 No    2+ races, Non-Hispanic                     25

ggplot(Q1.eth.table, aes(x=BMW1, y=count, fill=ppethm)) + geom_bar(position = "fill", stat="identity")

#household size
ggplot(Q1.comp, aes(x=BMW1, y=pphhsize)) + geom_boxplot()

#housing
Q1.hous.table=Q1.comp %>%
  group_by(BMW1, pphouse4) %>%
  summarize(count=n())
Q1.hous.table

## # A tibble: 8 × 3
## # Groups:   BMW1 [2]
##   BMW1  pphouse4                                              count
##   <fct> <fct>                                                 <int>
## 1 Yes   One-family house detached from any other house          306
## 2 Yes   One-family condo or townhouse attached to other units    19
## 3 Yes   Building with 2 or more apartments                       24
## 4 Yes   Other (mobile home, boat, RV, van, etc.)                 14
## 5 No    One-family house detached from any other house          414
## 6 No    One-family condo or townhouse attached to other units    72
## 7 No    Building with 2 or more apartments                      141
## 8 No    Other (mobile home, boat, RV, van, etc.)                 16

ggplot(Q1.hous.table, aes(x=BMW1, y=count, fill=pphouse4)) + geom_bar(position = "fill", stat="identity")

#income
Q1.inc.table=Q1.comp %>%
  group_by(BMW1, ppinc7) %>%
  summarize(count=n())
Q1.inc.table

## # A tibble: 14 × 3
## # Groups:   BMW1 [2]
##    BMW1  ppinc7               count
##    <fct> <fct>                <int>
##  1 Yes   Under $10,000            7
##  2 Yes   $10,000 to $24,999      24
##  3 Yes   $25,000 to $49,999      44
##  4 Yes   $50,000 to $74,999      65
##  5 Yes   $75,000 to $99,999      50
##  6 Yes   $100,000 to $149,999    80
##  7 Yes   $150,000 or more        93
##  8 No    Under $10,000           31
##  9 No    $10,000 to $24,999      53
## 10 No    $25,000 to $49,999     104
## 11 No    $50,000 to $74,999     105
## 12 No    $75,000 to $99,999      70
## 13 No    $100,000 to $149,999   115
## 14 No    $150,000 or more       165

ggplot(Q1.inc.table, aes(x=BMW1, y=count, fill=ppinc7)) + geom_bar(position = "fill", stat="identity")

#region
Q1.reg.table=Q1.comp %>%
  group_by(BMW1, ppreg4) %>%
  summarize(count=n())
Q1.reg.table

## # A tibble: 8 × 3
## # Groups:   BMW1 [2]
##   BMW1  ppreg4    count
##   <fct> <fct>     <int>
## 1 Yes   Northeast    43
## 2 Yes   Midwest      92
## 3 Yes   South       165
## 4 Yes   West         63
## 5 No    Northeast   139
## 6 No    Midwest     134
## 7 No    South       211
## 8 No    West        159

ggplot(Q1.reg.table, aes(x=BMW1, y=count, fill=ppreg4)) + geom_bar(position = "fill", stat="identity")

#kids under 17
Q1.kids=Q1.comp %>% 
  mutate(ppkid017 = case_when(ppkid017 <= 0 ~ "No",
                        ppkid017 > 0  ~ "Yes"))

Q1.kids.table=Q1.kids %>%
  group_by(BMW1, ppkid017) %>%
  summarize(count=n())
Q1.kids.table

## # A tibble: 4 × 3
## # Groups:   BMW1 [2]
##   BMW1  ppkid017 count
##   <fct> <chr>    <int>
## 1 Yes   No         252
## 2 Yes   Yes        111
## 3 No    No         469
## 4 No    Yes        174

Q1.kids.table$ppkid017=as.factor(Q1.kids.table$ppkid017)
ggplot(Q1.kids.table, aes(x=BMW1, y=count, fill=ppkid017)) + geom_bar(position = "fill", stat="identity")

#State
Q1.state.table=Q1.comp %>%
  group_by(BMW1, ppstaten) %>%
  summarize(count=n())
Q1.state.table=as.data.frame(Q1.state.table)
Q1.table.wide=reshape(Q1.state.table, idvar = "ppstaten",timevar = "BMW1", direction = "wide")
Q1.table.wide$Proportion=(Q1.table.wide$count.Yes/(Q1.table.wide$count.Yes+Q1.table.wide$count.No)*100)

Q1.table.wide=Q1.table.wide %>% rename(state = ppstaten)


plot_usmap(data=Q1.table.wide[,c(1,4)], regions = "states", values="Proportion") + 
  labs(title = "U.S. States",
       subtitle = "Proportion of gun owners by state.") + 
  theme(panel.background=element_blank())

BMW2. Safe Storage Views. (no to firearms)

#BMW2_1 BMW2_1. How much do you agree or disagree with each of the following statements?
  #People who have kids in the home should store their firearms safely so that the kids can’t access them.
#BMW2_2 BMW2_2. How much do you agree or disagree with each of the following statements?
  #Someone in the home struggling with mental health should not be able to access firearms.
#BMW2_3 BMW2_3. How much do you agree or disagree with each of the following statements?
  #Safe firearm storage is important to prevent my firearms from being stolen.
#BMW2_4 BMW2_4. How much do you agree or disagree with each of the following statements?
  #If I had a friend who I thought was considering suicide, I would ask them to lock up their firearms
#BMW2_5 BMW2_5. How much do you agree or disagree with each of the following statements?
  #Safe storage of firearms is important to prevent someone accidentally being hurt by them.
#BMW2_6 BMW2_6. How much do you agree or disagree with each of the following statements?
  #Safe storage of firearms is helpful for preventing them from being damaged (for example, in a fire).
#BMW2_7 BMW2_7. How much do you agree or disagree with each of the following statements?
  #It is important to keep at least one firearm unlocked and ready to use.
#BMW2_8 BMW2_8. How much do you agree or disagree with each of the following statements?
  #Quick access to firearms is an important way to keep my family safe.

itemsBMW2 <- data[, substr(names(data), 1, 5) == "BMW2_"] 
itemsBMW2[itemsBMW2 == "Skipped"] <- NA

names(itemsBMW2) <- c(
  BMW2_1="People who have kids in the home should store their firearms safely so that the kids can’t access them",
  BMW2_2="Someone in the home struggling with mental health should not be able to access firearms.",
  BMW2_3="Safe firearm storage is important to prevent my firearms from being stolen.",
  BMW2_4="If I had a friend who I thought was considering suicide, I would ask them to lock up their firearms",
  BMW2_5="Safe storage of firearms is important to prevent someone accidentally being hurt by them.",
  BMW2_6="Safe storage of firearms is helpful for preventing them from being damaged (for example, in a fire).",
  BMW2_7="It is important to keep at least one firearm unlocked and ready to use.",
  BMW2_8="Quick access to firearms is an important way to keep my family safe.")

p <- likert(itemsBMW2) 
plot(p)

plot(p, center=2, plot.percents=TRUE, plot.percent.low=FALSE, plot.percent.high=FALSE)

BMW2g <- likert(itemsBMW2, grouping=data$ppgender)
plot(BMW2g)

BMW2g2 <- likert(itemsBMW2, grouping=data$xparty4)
plot(BMW2g2)

BMW2g3 <- likert(itemsBMW2, grouping=data$xurbanicity)
plot(BMW2g3)

BMW2g4 <- likert(itemsBMW2, grouping=data$ppeducat)
plot(BMW2g4)

BMW2g5 <- likert(itemsBMW2, grouping=data$ppethm)
plot(BMW2g5)

BMW2g6 <- likert(itemsBMW2, grouping=data$ppinc7)
plot(BMW2g6)

BMW2g7 <- likert(itemsBMW2, grouping=data$ppmarit5)
plot(BMW2g7)

BMW2g8 <- likert(itemsBMW2, grouping=data$ppreg4)
plot(BMW2g8)

BMW3. How frequently do you carry a loaded firearm? (yes to firearms)

summary(data$BMW3)

##                                   Skipped 
##                                         7 
##                                     Daily 
##                                        37 
##   At least once a week, but not every day 
##                                        29 
## At least once a month, but not every week 
##                                        19 
## At least once a year, but not every month 
##                                        15 
##                              Almost never 
##                                        90 
## I own a firearm but never carry it loaded 
##                                       166 
##                                      NA's 
##                                       654

#State
Q3.comp=subset(data, BMW3 != "Skipped" | BMW3 != "NA")
summary(Q3.comp$BMW3)

##                                   Skipped 
##                                         7 
##                                     Daily 
##                                        37 
##   At least once a week, but not every day 
##                                        29 
## At least once a month, but not every week 
##                                        19 
## At least once a year, but not every month 
##                                        15 
##                              Almost never 
##                                        90 
## I own a firearm but never carry it loaded 
##                                       166

Q3.numer=data %>% 
    mutate(BMW3 = case_when(BMW3== "Daily" ~ "5.Daily", TRUE ~ BMW3))  %>%
   mutate(BMW3 = case_when(BMW3== "At least once a week, but not every day" ~ "4.Weekly", TRUE ~ BMW3))  %>%
   mutate(BMW3 = case_when(BMW3== "At least once a month, but not every week" ~ "3.Monthly", TRUE ~ BMW3))  %>%
   mutate(BMW3 = case_when(BMW3== "At least once a year, but not every month" ~ "2.Annually", TRUE ~ BMW3))  %>%
   mutate(BMW3 = case_when(BMW3== "Almost never" ~ "1.Almost never", TRUE ~ BMW3))  %>%
   mutate(BMW3 = case_when(BMW3== "I own a firearm but never carry it loaded" ~ "0.Never", TRUE ~ BMW3)) %>%
   mutate(BMW3 = case_when(BMW3== "Skipped" ~ "NA", TRUE ~ BMW3)) 

#Gender Subset

chisq.test(Q3.numer$BMW3, Q3.numer$ppgender)

## 
##  Pearson's Chi-squared test
## 
## data:  Q3.numer$BMW3 and Q3.numer$ppgender
## X-squared = 16.102, df = 6, p-value = 0.01322

Q3.table=Q3.numer %>%
  group_by(BMW3, ppgender) %>%
  filter(BMW3 != "NA") %>%
  summarize(count=n())
ggplot(Q3.table, aes(x=BMW3, y=count, fill=ppgender)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.table, aes(x=ppgender, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#party
Q3.party.table=Q3.numer %>%
  group_by(BMW3, xparty4) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.party.table

## # A tibble: 23 × 3
## # Groups:   BMW3 [6]
##    BMW3           xparty4        count
##    <chr>          <fct>          <int>
##  1 0.Never        Republican        58
##  2 0.Never        Democrat          38
##  3 0.Never        Independent       54
##  4 0.Never        Something else    16
##  5 1.Almost never Republican        36
##  6 1.Almost never Democrat          11
##  7 1.Almost never Independent       34
##  8 1.Almost never Something else     9
##  9 2.Annually     Republican         5
## 10 2.Annually     Democrat           1
## # ℹ 13 more rows

ggplot(Q3.party.table, aes(x=BMW3, y=count, fill=xparty4)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.party.table, aes(x=xparty4, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#urbanicity
Q3.urb.table=Q3.numer %>%
  group_by(BMW3, xurbanicity) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.urb.table

## # A tibble: 18 × 3
## # Groups:   BMW3 [6]
##    BMW3           xurbanicity count
##    <chr>          <fct>       <int>
##  1 0.Never        Urban          48
##  2 0.Never        Rural          40
##  3 0.Never        Suburban       78
##  4 1.Almost never Urban          25
##  5 1.Almost never Rural          26
##  6 1.Almost never Suburban       39
##  7 2.Annually     Urban           4
##  8 2.Annually     Rural           4
##  9 2.Annually     Suburban        7
## 10 3.Monthly      Urban           1
## 11 3.Monthly      Rural           5
## 12 3.Monthly      Suburban       13
## 13 4.Weekly       Urban           6
## 14 4.Weekly       Rural           6
## 15 4.Weekly       Suburban       17
## 16 5.Daily        Urban          11
## 17 5.Daily        Rural          15
## 18 5.Daily        Suburban       11

ggplot(Q3.urb.table, aes(x=BMW3, y=count, fill=xurbanicity)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.urb.table, aes(x=xurbanicity, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#age
ggplot(Q3.numer, aes(x=BMW3, y=ppage)) + geom_boxplot()

#education
Q3.ed.table=Q3.numer %>%
  group_by(BMW3, ppeducat) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.ed.table

## # A tibble: 23 × 3
## # Groups:   BMW3 [6]
##    BMW3           ppeducat           count
##    <chr>          <fct>              <int>
##  1 0.Never        Less than HS           9
##  2 0.Never        HS                    43
##  3 0.Never        Some college          48
##  4 0.Never        Bachelor or higher    66
##  5 1.Almost never Less than HS           3
##  6 1.Almost never HS                    24
##  7 1.Almost never Some college          33
##  8 1.Almost never Bachelor or higher    30
##  9 2.Annually     HS                     5
## 10 2.Annually     Some college           6
## # ℹ 13 more rows

ggplot(Q3.ed.table, aes(x=BMW3, y=count, fill=ppeducat)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.ed.table, aes(x=ppeducat, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#ethnicity
Q3.eth.table=Q3.numer %>%
  group_by(BMW3, ppethm) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.eth.table

## # A tibble: 24 × 3
## # Groups:   BMW3 [6]
##    BMW3           ppethm                                  count
##    <chr>          <fct>                                   <int>
##  1 0.Never        White, Non-Hispanic                       142
##  2 0.Never        Black or African American, Non-Hispanic     5
##  3 0.Never        Other, Non-Hispanic                         2
##  4 0.Never        Hispanic                                   12
##  5 0.Never        2+ races, Non-Hispanic                      5
##  6 1.Almost never White, Non-Hispanic                        67
##  7 1.Almost never Black or African American, Non-Hispanic     7
##  8 1.Almost never Other, Non-Hispanic                         3
##  9 1.Almost never Hispanic                                    8
## 10 1.Almost never 2+ races, Non-Hispanic                      5
## # ℹ 14 more rows

ggplot(Q3.eth.table, aes(x=BMW3, y=count, fill=ppethm)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.eth.table, aes(x=ppethm, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#household size
ggplot(Q3.numer, aes(x=BMW3, y=pphhsize)) + geom_boxplot()

#housing
Q3.hous.table=Q3.numer %>%
  group_by(BMW3, pphouse4) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.hous.table

## # A tibble: 20 × 3
## # Groups:   BMW3 [6]
##    BMW3           pphouse4                                              count
##    <chr>          <fct>                                                 <int>
##  1 0.Never        One-family house detached from any other house          140
##  2 0.Never        One-family condo or townhouse attached to other units    12
##  3 0.Never        Building with 2 or more apartments                       10
##  4 0.Never        Other (mobile home, boat, RV, van, etc.)                  4
##  5 1.Almost never One-family house detached from any other house           76
##  6 1.Almost never One-family condo or townhouse attached to other units     5
##  7 1.Almost never Building with 2 or more apartments                        5
##  8 1.Almost never Other (mobile home, boat, RV, van, etc.)                  4
##  9 2.Annually     One-family house detached from any other house           11
## 10 2.Annually     Building with 2 or more apartments                        3
## 11 2.Annually     Other (mobile home, boat, RV, van, etc.)                  1
## 12 3.Monthly      One-family house detached from any other house           16
## 13 3.Monthly      One-family condo or townhouse attached to other units     2
## 14 3.Monthly      Building with 2 or more apartments                        1
## 15 4.Weekly       One-family house detached from any other house           26
## 16 4.Weekly       Building with 2 or more apartments                        2
## 17 4.Weekly       Other (mobile home, boat, RV, van, etc.)                  1
## 18 5.Daily        One-family house detached from any other house           30
## 19 5.Daily        Building with 2 or more apartments                        3
## 20 5.Daily        Other (mobile home, boat, RV, van, etc.)                  4

ggplot(Q3.hous.table, aes(x=BMW3, y=count, fill=pphouse4)) + geom_bar(position = "fill", stat="identity")

#income
Q3.inc.table=Q3.numer %>%
  group_by(BMW3, ppinc7) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.inc.table

## # A tibble: 38 × 3
## # Groups:   BMW3 [6]
##    BMW3           ppinc7               count
##    <chr>          <fct>                <int>
##  1 0.Never        Under $10,000            2
##  2 0.Never        $10,000 to $24,999      16
##  3 0.Never        $25,000 to $49,999      16
##  4 0.Never        $50,000 to $74,999      30
##  5 0.Never        $75,000 to $99,999      21
##  6 0.Never        $100,000 to $149,999    32
##  7 0.Never        $150,000 or more        49
##  8 1.Almost never Under $10,000            1
##  9 1.Almost never $10,000 to $24,999       5
## 10 1.Almost never $25,000 to $49,999      13
## # ℹ 28 more rows

ggplot(Q3.inc.table, aes(x=BMW3, y=count, fill=ppinc7)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.inc.table, aes(x=ppinc7, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#region
Q3.reg.table=Q3.numer %>%
  group_by(BMW3, ppreg4) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.reg.table

## # A tibble: 23 × 3
## # Groups:   BMW3 [6]
##    BMW3           ppreg4    count
##    <chr>          <fct>     <int>
##  1 0.Never        Northeast    25
##  2 0.Never        Midwest      44
##  3 0.Never        South        65
##  4 0.Never        West         32
##  5 1.Almost never Northeast    11
##  6 1.Almost never Midwest      24
##  7 1.Almost never South        44
##  8 1.Almost never West         11
##  9 2.Annually     Midwest       1
## 10 2.Annually     South         9
## # ℹ 13 more rows

ggplot(Q3.reg.table, aes(x=BMW3, y=count, fill=ppreg4)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.reg.table, aes(x=ppreg4, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

#kids under 17
Q3.kids=Q3.numer %>% 
  mutate(ppkid017 = case_when(ppkid017 <= 0 ~ "No",
                        ppkid017 > 0  ~ "Yes"))

Q3.kids.table=Q3.kids %>%
  group_by(BMW3, ppkid017) %>%
    filter(BMW3 != "NA") %>%
  summarize(count=n())
Q3.kids.table

## # A tibble: 12 × 3
## # Groups:   BMW3 [6]
##    BMW3           ppkid017 count
##    <chr>          <chr>    <int>
##  1 0.Never        No         120
##  2 0.Never        Yes         46
##  3 1.Almost never No          72
##  4 1.Almost never Yes         18
##  5 2.Annually     No          12
##  6 2.Annually     Yes          3
##  7 3.Monthly      No           8
##  8 3.Monthly      Yes         11
##  9 4.Weekly       No          14
## 10 4.Weekly       Yes         15
## 11 5.Daily        No          23
## 12 5.Daily        Yes         14

Q3.kids.table$ppkid017=as.factor(Q3.kids.table$ppkid017)
ggplot(Q3.kids.table, aes(x=BMW3, y=count, fill=ppkid017)) + geom_bar(position = "fill", stat="identity")

ggplot(Q3.kids.table, aes(x=ppkid017, y=count, fill=BMW3)) + geom_bar(position = "fill", stat="identity")

(yes to firearms)

#BMW5   BMW5. Do you use any of the following measures to safely store your firearm(s)?
#columns 1,16-25,47-65

BMW5.only= data %>% 
  select(c(1, 16:25, 47:65)) %>%
  filter(BMW5_1 != "Skipped")

BMW5.only= BMW5.only %>% 
  mutate(BMW5_1 = case_when(BMW5_1== "Keeping all firearms unloaded" ~ "1", TRUE ~ BMW5_1))  %>%
  mutate(BMW5_2 = case_when(BMW5_2== "A safe, lockbox, or gun case (with a key, access code, or biometric lock)" ~ "1", TRUE ~ BMW5_2))  %>%
  mutate(BMW5_3 = case_when(BMW5_3== "A cable or trigger lock" ~ "1", TRUE ~ BMW5_3))  %>%
  mutate(BMW5_4 = case_when(BMW5_4== "Motion alarm or notifications for firearm access" ~ "1", TRUE ~ BMW5_4))  %>%
  mutate(BMW5_5 = case_when(BMW5_5== "Firearm sensor for close range or self-aiming" ~ "1", TRUE ~ BMW5_5))  %>%
  mutate(BMW5_6 = case_when(BMW5_6== "Locking up ammunition separate from firearms" ~ "1", TRUE ~ BMW5_6))  %>%
  mutate(BMW5_7 = case_when(BMW5_7== "Disassembling firearms and storing parts separately" ~ "1", TRUE ~ BMW5_7))  %>%
  mutate(BMW5_8 = case_when(BMW5_8== "Entrusting keys or firearms to a trusted person or organization" ~ "1", TRUE ~ BMW5_8))  %>%
  mutate(BMW5_9 = case_when(BMW5_9== "Using other safety measure not listed here (please specify)" ~ "1", TRUE ~ BMW5_9))  %>%
  mutate(BMW5_10 = case_when(BMW5_10== "None of the above" ~ "1", TRUE ~ BMW5_10))

BMW5.only <- BMW5.only %>% mutate_at(c(2:11), as.numeric)

str(BMW5.only)

## 'data.frame':    355 obs. of  30 variables:
##  $ CaseID     : num  3 8 10 11 13 14 19 24 25 26 ...
##  $ BMW5_1     : num  0 1 0 1 0 1 0 0 0 0 ...
##  $ BMW5_2     : num  0 1 1 0 0 1 0 0 1 0 ...
##  $ BMW5_3     : num  0 1 1 0 0 0 0 0 0 0 ...
##  $ BMW5_4     : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ BMW5_5     : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ BMW5_6     : num  0 1 0 0 0 1 0 1 1 0 ...
##  $ BMW5_7     : num  0 0 0 1 0 0 0 0 0 0 ...
##  $ BMW5_8     : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ BMW5_9     : num  0 0 0 1 0 0 0 0 0 0 ...
##  $ BMW5_10    : num  1 0 0 0 1 0 1 0 0 1 ...
##  $ xparty4    : Factor w/ 5 levels "Republican","Democrat",..: 1 3 3 2 3 1 1 3 1 3 ...
##  $ xurbanicity: Factor w/ 3 levels "Urban","Rural",..: 1 3 3 3 1 1 3 1 2 3 ...
##  $ ppage      : num  65 56 75 66 79 23 56 56 71 75 ...
##  $ ppeducat   : Factor w/ 4 levels "Less than HS",..: 3 4 4 4 2 3 4 3 3 4 ...
##  $ ppeduc5    : Factor w/ 5 levels "No high school diploma or GED",..: 3 5 4 4 2 3 4 3 3 5 ...
##  $ ppethm     : Factor w/ 5 levels "White, Non-Hispanic",..: 1 1 1 1 1 1 1 1 1 2 ...
##  $ ppgender   : Factor w/ 2 levels "Male","Female": 1 1 1 1 2 1 2 2 2 1 ...
##  $ pphhsize   : num  2 4 7 3 1 2 4 2 2 2 ...
##  $ pphouse4   : Factor w/ 4 levels "One-family house detached from any other house",..: 1 1 1 1 1 4 1 1 1 1 ...
##  $ ppinc7     : Factor w/ 7 levels "Under $10,000",..: 4 7 5 6 5 4 5 3 6 7 ...
##  $ ppmarit5   : Factor w/ 5 levels "Now married",..: 1 1 1 1 2 5 3 2 1 1 ...
##  $ ppmsacat   : Factor w/ 2 levels "Non-metro area",..: 2 2 2 2 2 2 2 2 2 2 ...
##  $ ppreg4     : Factor w/ 4 levels "Northeast","Midwest",..: 3 3 3 3 3 4 3 3 2 2 ...
##  $ pprent     : Factor w/ 3 levels "Owned or being bought by you or someone in",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ ppemploy   : Factor w/ 3 levels "Working full-time",..: 3 1 3 3 3 3 1 1 2 3 ...
##  $ ppstaten   : Factor w/ 51 levels "ME","NH","VT",..: 30 28 30 28 38 47 30 38 14 12 ...
##  $ ppt18ov    : num  2 3 6 2 1 2 4 2 2 2 ...
##  $ ppkid017   : num  0 1 1 1 0 0 0 0 0 0 ...
##  $ PPREG9     : Factor w/ 9 levels "New England",..: 5 5 5 5 7 9 5 7 3 3 ...
##  - attr(*, "variable.labels")= Named chr [1:67] "CaseID" "Response ID" "Start Date" "End Date" ...
##   ..- attr(*, "names")= chr [1:67] "CaseID" "ResponseId" "StartDate" "EndDate" ...
##  - attr(*, "codepage")= int 1252

BMW5.only=BMW5.only %>% 
  rename(
    unloaded = BMW5_1,
    safe_case_lockbox = BMW5_2,
    cable_trigger.lock = BMW5_3,
    access.alarm = BMW5_4,
    range.sensor = BMW5_5,
    sep.ammunition = BMW5_6,
    disassembled_sep.storage = BMW5_7,
    trusted.person_org = BMW5_8,
    other = BMW5_9,
    none.of.above = BMW5_10
    )



BMW5.perc=((colSums(BMW5.only[,2:11], na.rm = T))/355)*100
BMW5.perc=as.data.frame(BMW5.perc, row.names = F)

BMW5.perc <- rownames_to_column(BMW5.perc, var = "Method")
names(BMW5.perc)[2] <- "Percentage"

BMW5.perc

##                      Method Percentage
## 1                  unloaded  51.267606
## 2         safe_case_lockbox  59.154930
## 3        cable_trigger.lock  21.690141
## 4              access.alarm   2.816901
## 5              range.sensor   1.126761
## 6            sep.ammunition  38.309859
## 7  disassembled_sep.storage   6.760563
## 8        trusted.person_org   0.000000
## 9                     other   3.098592
## 10            none.of.above  14.929577

ggplot(BMW5.perc, aes(reorder(x=Method , -Percentage), y=Percentage)) + geom_bar( stat="identity") +
  labs(title="Current Safety Precautions Used")

    #1 Keeping all firearms unloaded
    #2 A safe, lockbox, or gun case (with a key, access code, or biometric lock)
    #3 A cable or trigger lock
    #4 Motion alarm or notifications for firearm access
    #5 Firearm sensor for close range or self-aiming
    #6 Locking up ammunition separate from firearms
    #7 Disassembling firearms and storing parts separately
    #8 Entrusting keys or firearm parts to a trusted person or organization
    #9 Using other safety measure not listed here
    #10 None of the above

(everyone)

#BMW6   BMW6. Which of the following methods for storing firearms that you are not currently using/would you consider using?

BMW6.only= data %>% 
  select(c(1, 27:36, 47:65)) %>%
  filter(BMW6_1 != "Skipped")

BMW6.only= BMW6.only %>% 
  mutate(BMW6_1 = case_when(BMW6_1== "Keeping all firearms unloaded" ~ "1", TRUE ~ BMW6_1))  %>%
  mutate(BMW6_2 = case_when(BMW6_2== "A safe, lockbox, or gun case (with a key, access code, or biometric lock)" ~ "1", TRUE ~ BMW6_2))  %>%
  mutate(BMW6_3 = case_when(BMW6_3== "A cable or trigger lock" ~ "1", TRUE ~ BMW6_3))  %>%
  mutate(BMW6_4 = case_when(BMW6_4== "Motion alarm or notifications for firearm access" ~ "1", TRUE ~ BMW6_4))  %>%
  mutate(BMW6_5 = case_when(BMW6_5== "Firearm sensor for close range or self-aiming" ~ "1", TRUE ~ BMW6_5))  %>%
  mutate(BMW6_6 = case_when(BMW6_6== "Locking up ammunition separate from firearms" ~ "1", TRUE ~ BMW6_6))  %>%
  mutate(BMW6_7 = case_when(BMW6_7== "Disassembling firearms and storing parts separately" ~ "1", TRUE ~ BMW6_7))  %>%
  mutate(BMW6_8 = case_when(BMW6_8== "Entrusting keys or firearms to a trusted person or organization" ~ "1", TRUE ~ BMW6_8))  %>%
  mutate(BMW6_9 = case_when(BMW6_9== "Using other safety measure not listed here (please specify)" ~ "1", TRUE ~ BMW6_9))  %>%
  mutate(BMW6_10 = case_when(BMW6_10== "None of the above" ~ "1", TRUE ~ BMW6_10))

BMW6.only <- BMW6.only %>% mutate_at(c(2:11), as.numeric)

str(BMW6.only)

## 'data.frame':    811 obs. of  30 variables:
##  $ CaseID     : num  1 2 3 4 5 6 7 9 10 12 ...
##  $ BMW6_1     : num  0 0 1 1 1 0 1 0 1 1 ...
##  $ BMW6_2     : num  1 1 1 1 1 1 1 1 NA 1 ...
##  $ BMW6_3     : num  1 0 0 1 0 1 1 1 NA 0 ...
##  $ BMW6_4     : num  0 0 0 1 0 0 0 0 0 0 ...
##  $ BMW6_5     : num  0 0 0 0 0 1 0 0 0 0 ...
##  $ BMW6_6     : num  NA 0 0 NA NA NA NA NA NA NA ...
##  $ BMW6_7     : num  0 0 0 0 1 0 1 0 1 0 ...
##  $ BMW6_8     : num  0 0 0 0 0 0 0 0 NA 0 ...
##  $ BMW6_9     : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ BMW6_10    : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ xparty4    : Factor w/ 5 levels "Republican","Democrat",..: 3 4 1 3 3 2 2 3 3 1 ...
##  $ xurbanicity: Factor w/ 3 levels "Urban","Rural",..: 2 1 1 1 1 2 1 1 3 2 ...
##  $ ppage      : num  64 66 65 61 41 71 60 70 75 75 ...
##  $ ppeducat   : Factor w/ 4 levels "Less than HS",..: 1 3 3 4 3 4 4 4 4 3 ...
##  $ ppeduc5    : Factor w/ 5 levels "No high school diploma or GED",..: 1 3 3 4 3 5 4 4 4 3 ...
##  $ ppethm     : Factor w/ 5 levels "White, Non-Hispanic",..: 5 1 1 1 3 1 3 1 1 1 ...
##  $ ppgender   : Factor w/ 2 levels "Male","Female": 2 1 1 2 1 2 2 1 1 2 ...
##  $ pphhsize   : num  2 1 2 3 4 2 2 2 7 2 ...
##  $ pphouse4   : Factor w/ 4 levels "One-family house detached from any other house",..: 1 1 1 2 3 1 3 1 1 1 ...
##  $ ppinc7     : Factor w/ 7 levels "Under $10,000",..: 7 6 4 7 3 7 7 7 5 3 ...
##  $ ppmarit5   : Factor w/ 5 levels "Now married",..: 1 2 1 1 1 1 1 1 1 1 ...
##  $ ppmsacat   : Factor w/ 2 levels "Non-metro area",..: 2 2 2 2 2 1 2 1 2 1 ...
##  $ ppreg4     : Factor w/ 4 levels "Northeast","Midwest",..: 3 4 3 2 1 1 1 3 3 3 ...
##  $ pprent     : Factor w/ 3 levels "Owned or being bought by you or someone in",..: 1 2 1 1 2 1 1 1 1 1 ...
##  $ ppemploy   : Factor w/ 3 levels "Working full-time",..: 3 3 3 1 1 3 1 1 3 3 ...
##  $ ppstaten   : Factor w/ 51 levels "ME","NH","VT",..: 35 49 30 10 7 6 7 27 30 31 ...
##  $ ppt18ov    : num  2 1 2 3 2 2 2 2 6 2 ...
##  $ ppkid017   : num  0 0 0 0 2 0 0 0 1 0 ...
##  $ PPREG9     : Factor w/ 9 levels "New England",..: 7 9 5 3 2 1 2 5 5 6 ...
##  - attr(*, "variable.labels")= Named chr [1:67] "CaseID" "Response ID" "Start Date" "End Date" ...
##   ..- attr(*, "names")= chr [1:67] "CaseID" "ResponseId" "StartDate" "EndDate" ...
##  - attr(*, "codepage")= int 1252

BMW6.only=BMW6.only %>% 
  rename(
    unloaded = BMW6_1,
    safe_case_lockbox = BMW6_2,
    cable_trigger.lock = BMW6_3,
    access.alarm = BMW6_4,
    range.sensor = BMW6_5,
    sep.ammunition = BMW6_6,
    disassembled_sep.storage = BMW6_7,
    trusted.person_org = BMW6_8,
    other = BMW6_9,
    none.of.above = BMW6_10
    )

BMW6.perc=((colSums(BMW6.only[,2:11], na.rm = T))/811)*100
BMW6.perc=as.data.frame(BMW6.perc, row.names = F)

BMW6.perc <- rownames_to_column(BMW6.perc, var = "Method")
names(BMW6.perc)[2] <- "Percentage"

BMW6.perc

##                      Method Percentage
## 1                  unloaded  47.348952
## 2         safe_case_lockbox  68.803946
## 3        cable_trigger.lock  34.278668
## 4              access.alarm  28.113440
## 5              range.sensor  14.056720
## 6            sep.ammunition   0.000000
## 7  disassembled_sep.storage  16.892725
## 8        trusted.person_org   0.000000
## 9                     other   3.822441
## 10            none.of.above  16.892725

ggplot(BMW6.perc, aes(reorder(x=Method , -Percentage), y=Percentage)) + geom_bar( stat="identity") +
  labs(title="Current Safety Precautions Used")

   #1 Keeping all firearms unloaded
   #2 A safe, lockbox, or gun case (with a key, access code, or biometric lock)
   #3 A cable or trigger lock
   #4  Motion alarm or notifications for firearm access
   #5 Firearm sensor for close range or self-aiming
   #6 Locking up ammunition separate from firearms
   #7 Disassembling firearms and storing parts separately
   #8 Entrusting keys or firearm parts to a trusted person or organization
   #9 Using other safety measure not listed here
   #10 None of the above
   #9_text

names(BMW5.perc)[2] <- "Percentage.Used"
names(BMW6.perc)[2] <- "Percentage.Consider"

merged_methods <- merge(BMW5.perc, BMW6.perc, by = "Method")

long_data <- reshape(merged_methods, 
                     varying = list(c("Percentage.Used", "Percentage.Consider")), 
                     direction = "long", 
                     times = c("Percentage.Used", "Percentage.Consider"), 
                     v.names = "Percentage", 
                     timevar = "Status", 
                     idvar = "Method")

ggplot(long_data, aes(fill=Status, y=Percentage, reorder(x=Method , -Percentage))) + 
    geom_bar(position="dodge", stat="identity")

(yes to firearms)

#BMW4_1 BMW4_1. How much do you agree or disagree with each of the following statements?
  #If I had kids at home I would lock up my firearms
#BMW4_2 BMW4_2. How much do you agree or disagree with each of the following statements?
  #If someone in my home was struggling with mental or physical health, I would lock up my firearm(s).
#BMW4_3 BMW4_3. How much do you agree or disagree with each of the following statements?
  #If a close friend or family asked me to, I would lock up my firearm(s).
#BMW4_4 BMW4_4. How much do you agree or disagree with each of the following statements?
  #To prevent my firearm(s) from being used, stolen, or damaged, I would lock it up.
#BMW4_5 BMW4_5. How much do you agree or disagree with each of the following statements?
  #To prevent accidental injury from my firearm(s), I would lock it up.
#BMW4_6 BMW4_6. How much do you agree or disagree with each of the following statements?
  #To prevent suicide for me or others, I would lock up my firearm(s).
#BMW4_7 BMW4_7. How much do you agree or disagree with each of the following statements?
  #I plan to always keep at least one firearm unlocked.
itemsBMW4 <- data[, substr(names(data), 1, 5) == "BMW4_"] 
itemsBMW4[itemsBMW4 == "Skipped"] <- NA

names(itemsBMW4) <- c(
  BMW4_1="If I had kids at home I would lock up my firearms",
  BMW4_2="If someone in my home was struggling with mental or physical health, I would lock up my firearm(s).",
  BMW4_3="If a close friend or family asked me to, I would lock up my firearm(s).",
  BMW4_4="To prevent my firearm(s) from being used, stolen, or damaged, I would lock it up.",
  BMW4_5="To prevent accidental injury from my firearm(s), I would lock it up.",
  BMW4_6="To prevent suicide for me or others, I would lock up my firearm(s).",
  BMW4_7="I plan to always keep at least one firearm unlocked.")

p <- likert(itemsBMW4) 
plot(p)

plot(p, center=2, plot.percents=TRUE, plot.percent.low=FALSE, plot.percent.high=FALSE)

BMW4g <- likert(itemsBMW4, grouping=data$ppgender)
plot(BMW4g)

BMW4g2 <- likert(itemsBMW4, grouping=data$xparty4)
plot(BMW4g2)

BMW4g3 <- likert(itemsBMW4, grouping=data$xurbanicity)
plot(BMW4g3)

BMW4g4 <- likert(itemsBMW4, grouping=data$ppeducat)
plot(BMW4g4)

BMW4g5 <- likert(itemsBMW4, grouping=data$ppethm)
plot(BMW4g5)

BMW4g6 <- likert(itemsBMW4, grouping=data$ppinc7)
plot(BMW4g6)

BMW4g7 <- likert(itemsBMW4, grouping=data$ppmarit5)
plot(BMW4g7)

BMW4g8 <- likert(itemsBMW4, grouping=data$ppreg4)
plot(BMW4g8)

BMW4g9 <- likert(itemsBMW4, grouping=data$ppkid017)
plot(BMW4g9)

Quick access to firearm #2 (all responses) marital status by gender

hist(Q3.numer$BMW3) # See number of individuals of male/female in each category

Black people and their daily carrying (1/3 carrying weekly) - same thing for income - look at sample size for lowest income - look at income by race category

-crosstab between region and urban/suburban/rural

Q3.party.table

people with kids are more commonly carrying? Q3

For comparison of 5/6 restrict to only people who answered both