Bootstrapping

Learning Objectives

At the end of this lab, you will:

1. Understand the principles of bootstrapping
2. Understand how to apply the bootstrap confidence interval to inference in linear models

What You Need

1. Be up to date with lectures
2. Have completed previous lab exercises from Week 2 and Week 9

Required R Packages

Remember to load all packages within a code chunk at the start of your RMarkdown file using library(). If you do not have a package and need to install, do so within the console using install.packages(" "). For further guidance on installing/updating packages, see Section C here.

For this lab, you will need to load the following package(s):

• tidyverse
• psych
• patchwork
• sjPlot
• kableExtra
• car

Presenting Results

All results should be presented following APA guidelines. If you need a reminder on how to hide code, format tables/plots, etc., make sure to review the rmd bootcamp.

The example write-up sections included as part of the solutions are not perfect - they instead should give you a good example of what information you should include and how to structure this. Note that you must not copy any of the write-ups included below for future reports - if you do, you will be committing plagiarism, and this type of academic misconduct is taken very seriously by the University. You can find out more here.

Lab Data

You can download the data required for this lab here or read it in via this link https://uoepsy.github.io/data/science-faith-attitude.csv

Study Overview

Research Question

Is there an association between peoples’ attitudes towards science and faith and their scientific knowledge after accounting for their age?

Attitudes Codebook

Description

Using Eurobarometer measures, researchers were interested in understanding opinions about science, technology, social values, and services of general interest. Specifically, researchers wanted to explore whether peoples’ attitudes towards science and faith was associated with their scientific knowledge after accounting for their age. To examine these associations, we will use a subset of data from the 2005 Eurobarometer 63.1 survey.

Data Dictionary

variable	description
v5	Participant ID
v6	Nation ID: 1 = France, 2 = Belgium, 3 = The Netherlands, 4 = Germany West, 5 = Italy, 6 = Luxembourg, 7 = Denmark, 8 = Ireland, 9 = Great Britain, 10 = Northern Ireland, 11 = Greece, 12 = Spain, 13 = Portugal, 14 = Germany East, 15 = Norway, 16 = Finland, 17 = Sweden, 18 = Austria, 19 = Cyprus (Republic), 20 = Czech Republic, 21 = Estonia, 22 = Hungary, 23 = Latvia, 24 = Lithuania, 25 = Malta, 26 = Poland, 27 = Slovakia, 28 = Slovenia, 29 = Bulgaria, 30 = Romania, 31 = Turkey, 32 = Croatia, 33 Cyprus (TCC; not included), 34 = Iceland, 35 = Switzerland
kstot	Score on a science quiz composed of 13 true/false items. Scores can range 0-13.
age	Individual's age (in years)
male	Whether participant idenitfied as male (1) or female (0)
toomuchscience	Attitude to science and faith. Response to the question - We rely too much on science and not enough on faith. Responses were recorded on a 5-point scale from strongly disagree (0) to strongly agree (4)
solveprob	Number of explicitly mentioned target phrases in response to the question - Please tell me, in your own words, what it means to study something scientifically?

Preview

The first six rows of the data are:

v5	v6	kstot	age	male	toomuchscience	solveprob
1	2	12	35	1	2	0
2	2	11	34	0	1	1
3	2	7	40	0	NA	NA
4	2	9	32	0	3	0
5	2	6	35	0	3	1
6	2	11	67	1	NA	NA

Setup

Setup

1. Create a new RMarkdown file
2. Load the required package(s)
3. Read in the science-faith-attitude dataset into R, assigning it to an object named ebsurvey

Solution

#load packages
library(tidyverse) 
library(psych)
library(patchwork)
library(sjPlot)
library(kableExtra)
library(car)

#read in data
ebsurvey <- read_csv("https://uoepsy.github.io/data/science-faith-attitude.csv")

Exercises

Study Overview & Data Management

Question 1

Examine the dataset, and perform any necessary and appropriate data management steps.

Note, to address the research question, we only need to refer to the kstot, age, and toomuchscience variables. Subset the data to only have those 3 columns.

Hint

• Check for missing data. There are numerous ways to deal with this. Two common commands are na.omit() and drop_na(). The former will remove all rows from a dataset that contain NA values in any column. In the latter, we can specify which columns we want to identify NA values in, and remove only rows containing NA values for those specific columns. In other words, the latter can help to preserve more data.
• If needed, provide better variable names
• To subset the data to only include the 3 variables of interest, we can use the select() function

Solution

# Inspect top 6 rows
head(ebsurvey)

# A tibble: 6 × 7
     v5    v6 kstot   age  male toomuchscience solveprob
  <dbl> <dbl> <dbl> <dbl> <dbl>          <dbl>     <dbl>
1     1     2    12    35     1              2         0
2     2     2    11    34     0              1         1
3     3     2     7    40     0             NA        NA
4     4     2     9    32     0              3         0
5     5     2     6    35     0              3         1
6     6     2    11    67     1             NA        NA

# Check data dimensions
dim(ebsurvey)

[1] 21886     7

There are 21886 observations on 7 variables.

However, today we will be only using the kstot, age, and toomuchscience variables, and so we subset the data to only include these:

ebsurvey <- ebsurvey %>%
    select(kstot, age, toomuchscience)

Are there any NA values in the data?

#check for NAs
anyNA(ebsurvey)

[1] TRUE

#how many NAs?
table(is.na(ebsurvey))

FALSE  TRUE 
54268 11390 

#11390 NAs in data set 

#Omit the NAs - we are interested in all three columns so do not need to specify within drop_na()
ebsurvey <- ebsurvey %>%
    drop_na()

# Check new data dimensions
dim(ebsurvey)

[1] 10503     3

Give the variables more meaningful names. Rename kstot to science_knowledge and rename toomuchscience to attitude:

ebsurvey <- ebsurvey %>%
    rename(science_knowledge = kstot,
           attitude = toomuchscience)
head(ebsurvey)

# A tibble: 6 × 3
  science_knowledge   age attitude
              <dbl> <dbl>    <dbl>
1                12    35        2
2                11    34        1
3                 9    32        3
4                 6    35        3
5                 9    37        1
6                 5    63        2

Question 2

Provide a brief overview of the study design and data, before detailing your analysis plan to address the research question.

Hint

• Give the reader some background on the context of the study
• State what type of analysis you will conduct in order to address the research question
• Specify the model to be fitted to address the research question
• Specify your chosen significance (\(\alpha\)) level
• State your hypotheses

Much of the information required can be found in the Study Overview codebook.

Solution

The ebsurvey dataset, excluding missing values, included 10,503 individual respondents who were measured on 3 different attributes of interest: (1) score on a science “quiz” composed of 13 true/false items; (2) attitudes towards science and faith (question phrasing: “We rely too much on science and not enough on faith” with responses recorded on a 5-point scale from strongly disagree to strongly agree); and (3) their age (in years).

Histograms will be used to visualise the marginal distributions of attitudes towards science and faith, scientific knowledge, and age. To understand the strength of association among the variables, we will estimate the the correlation coefficients. To address the research question of whether there is an association between people’s scientific knowledge and their attitudes towards science and faith after accounting for their age, we are going to fit the following multiple linear regression model:

\[ \text{Attitude} = \beta_0 + \beta_1 \cdot \text{Science Knowledge} + \beta_2 \cdot \text{Age} + \epsilon \] Effects will be considered statistically significant at \(\alpha = .05\).

Our hypotheses are:

\(H_0: \beta_1 = 0\): There is no association between people’s scientific knowledge and their attitudes towards science and faith after accounting for their age

\(H_1: \beta_1 \neq 0\): There is an association between people’s scientific knowledge and their attitudes towards science and faith after accounting for their age

Descriptive Statistics & Visualisations

Question 3

Alongside descriptive statistics, visualize the marginal distributions of the attitude, science_knowledge, and age variables.

Hint

Code for tables & plots
- You might be able to re-use some of the code you wrote for Lab 1 here, and remember that you can refer back to the DAPR1 materials too

Plot interpretation
- The shape, center and spread of the distribution
- Whether the distribution is symmetric or skewed
- Whether the distribution is unimodal or bimodal

Plotting tips
- Use \n to wrap text in your titles and or axis labels
- You may wish to consider using helpful bin sizes, see the binwidth argument, but this is subjective and optional - The patchwork package allows us to arrange multiple plots in two ways - | arranges the plots adjacent to one another, and / arranges the plots on top of one another

Table tips
- The describe() function from the psych package will produce a table of descriptive statistics. If you would like only a subset of this output (e.g., mean, sd), you can use select() after calling describe() e.g., describe() %>% select(mean, sd)
- The kableExtra package allows us to produce well formatted tables for our descriptive statistics. To do so, you need to specify the kable() and kable_styling() arguments
- Review the guidance on the rmd bootcamp, particularly Lesson 4

Solution

Numeric

# drop columns 5, 6, 7 and 13 from describe output
ebsurvey %>% 
    describe() %>%
    select(-c(5:7, 13)) %>%
    kable(caption = "Attitude, Scientfic Knowledge, and Age Descriptive Statistics", align = "c", digits = 2) %>%
    kable_styling(full_width = FALSE)

Table 1: Attitude, Scientfic Knowledge, and Age Descriptive Statistics

	vars	n	mean	sd	min	max	range	skew	kurtosis
science_knowledge	1	10503	8.68	2.61	0	13	13	-0.42	-0.38
age	2	10503	44.93	17.32	15	93	78	0.18	-0.88
attitude	3	10503	2.20	1.22	0	4	4	-0.25	-0.87

Visual

Attitude Towards Science & Faith

ggplot(ebsurvey, aes(x = attitude)) +
    geom_bar() +
    labs(x = 'We Rely too Much on Science and not Enough on Faith', 
         y = 'Frequency')

Figure 1: Bar Chart of Attitude Scores

The mean score on the science and faith attitude variable is just over 2. There are only 5 discrete values possible in the distribution, based on the response options available, but the distribution looks approximately normal, with a slight negative skew (see Figure 1).

Science Knowledge

ggplot(ebsurvey, aes(x = science_knowledge)) +
    geom_histogram() +
    labs(x = 'Science Knowledge Quiz Scores', 
         y = 'Frequency')

Figure 2: Histogram of Scientific Knowledge Scores

Figure 2 shows that the majority of values on the science knowledge quiz score cluster between about 5 and 11. There is a slight negative skew to the distribution. Overall there is little reason for concern as to the appropriateness of the variable for inclusion.

Age

ggplot(ebsurvey, aes(x = age)) +
    geom_histogram() +
    labs(x = 'Age (years)', 
         y = 'Frequency')

Figure 3: Histogram of Age

The mean age in the sample is about 45 years with a standard deviation of just over 17 years. The distribution looks approximately normal, with a slight positive skew (see Figure 3).

Question 4

Produce plots of the associations between the outcome variable and each of the explanatory variables.

Hint

Using geom_point() here might not be the best idea, all we would see from the plot is all combinations of Attitude towards Science and Faith and a) Age and b) Science Knowledge that were observed in the data, and this is not very informative:

Figure 4: Scatterplots displaying the associations between Attitude towards Science and Faith and a) Age, and b) Science Knowledge

Instead, you may want to consider using geom_jitter() to add a little bit of noise (or jitter) to the plot. Within the geom_jitter() argument, take some time to experiment with the size = and alpha = arguments to find optimal values to aid interpretation.

Solution

p1 <- ggplot(data = ebsurvey, aes(x = age, y = attitude)) +
    geom_jitter(size = .5, alpha = .3) +
    labs(x = 'Age (in years)', y = "Attitude towards Science and Faith")

p2 <- ggplot(data = ebsurvey, aes(x = science_knowledge, y = attitude)) +
    geom_jitter(size = .5, alpha = .3)  +
    labs(x = "Science Knowledge Quiz Scores", y = "Attitude towards Science and Faith")

p1 | p2 

Figure 5: Scatterplots displaying the associations between Attitude towards Science and Faith and a) Age, and b) Science Knowledge

It does not seem like there is a strong linear dependence of attitude to science and faith on a person’s age or their scientific knowledge. We can make some rough observations:

• The majority of respondents scored 2 or 3 (i.e., responded as neutral or agree) on the attitudes towards science and faith question
• Very few people strongly disagreed with the attitudes towards science and faith question, but this appeared to be slightly more common in younger respondents
• Very few people had scientific knowledge quiz scores <5. Lower scientific knowledge quiz scores appeared to be associated with responses of 3 or 4 (i.e., responded as agree or strongly agree) on the attitudes towards science and faith question.

Question 5

Produce a correlation matrix of the variables which are to be used in the analysis, and write a short paragraph describing the associations.

Hint

Correlation Matrix
Review Lab 1 Q2 for guidance on how to produce a correlation matrix.

Interpret Strength of Correlation See DAPR1 Semester 2 Week 10 Lectures for a recap on how to interpret \(r\) values.

APA Format
Make sure to round your numbers in-line with APA 7th edition guidelines as noted at the start of the lab (see ‘Presenting Results’). The round() function will come in handy here.

Solution

We can either index the dataframe or select the variables of interest:

Index dataframe ([])

# correlation matrix of the three columns of interest (check which columns we need - in this case, 1,2,3)
round(cor(ebsurvey[,c(1:3)]), digits = 2)

                  science_knowledge   age attitude
science_knowledge              1.00 -0.12    -0.18
age                           -0.12  1.00     0.05
attitude                      -0.18  0.05     1.00

Variable selection (select())

# select only the columns we want by variable name, and pass this to cor()
ebsurvey %>% 
  select(attitude, science_knowledge, age) %>%
  cor() %>%
    kable(digits = 2, caption = "Correlation Matrix") %>%
    kable_styling(full_width = FALSE)

Table 2: Correlation Matrix

	attitude	science_knowledge	age
attitude	1.00	-0.18	0.05
science_knowledge	-0.18	1.00	-0.12
age	0.05	-0.12	1.00

• There was a weak, positive, linear association between attitude towards science and faith and age for the participants in the sample \((r = .05)\)
• There was a weak, negative, linear association between attitude towards science and faith and scientific knowledge for the participants in the sample \((r = -.18)\)
• There was a weak, negative, linear association between scientific knowledge and age for the participants in the sample \((r = -.12)\). The correlation is relatively small in absolute terms, and we therefore have little concern about multicollinearity influencing this regression analysis
• Overall, there were very weak linear associations among the variables of interest

Model Fitting & Interpretation

Question 6

Fit the specified model, and assign it the name “att_mdl”.

\[ \text{Attitude} = \beta_0 + \beta_1 \cdot \text{Science Knowledge} + \beta_2 \cdot \text{Age} + \epsilon \]

Solution

att_mdl <- lm(attitude ~ science_knowledge + age, data = ebsurvey)

Question 7

Check the assumptions of your model. Note any violations of the model assumptions.

Hint

The most efficient way to do this is likely to use the plot(model) function.

To have all 4 plots display in a nice 2 by 2 matrix, we can use the command par(mfrow = c(2,2)). Then create the plot. Make sure to always re-set this back to a single figure made up by a single panel with the command par(mfrow = c(1,1)).

Solution

par(mfrow = c(2,2)) # set 2 by 2 panels
plot(att_mdl)
par(mfrow = c(1,1)) # go back to 1 by 1 panels

Figure 6: Assumption Plots for Attitudes Model

Based on the visual inspection of the plots, the assumptions appear to be violated.

Question 8

Bootstrap your model, computing 1000 bootstrap samples.

Provide key model results in a formatted table, and interpret your coefficients in the context of the research question.

Hint

You can’t use tab_model() here, but instead will need to use kable() and kable_styling().

Solution

#Run model
boot_mdl <- Boot(att_mdl, R = 1000)

summary(boot_mdl) %>%
    kable(digits = 3, caption = 'Bootstrap Regression Results for Attitudes Model') %>%
    kable_styling(full_width = FALSE)

Table 3: Bootstrap Regression Results for Attitudes Model

	R	original	bootBias	bootSE	bootMed
(Intercept)	1000	2.788	0.002	0.053	2.790
science_knowledge	1000	-0.080	0.000	0.004	-0.080
age	1000	0.002	0.000	0.001	0.002

Intercept

Hint - Interpretation

The intercept of a multiple regression model can be interpreted as the average expected value of the dependent variable when all of the independent variables equal zero.

• As presented in Table 3, the estimated intercept was approximately 2.79. This represents the expected attitude towards science and faith score when all other variables were zero.
• In this case, only a handful of respondents had a scientific knowledge quiz score of zero, and nobody was aged zero. Hence, in this example, the intercept itself isn’t very informative.
  • A way to improve the model here to aid interpretation would be to mean centre age, and then refit the model with the mean centred age variable!

Scientific knowledge

Hint - Interpretation

This represents the average marginal effect of \(X1\) on \(Y\), and can be interpreted as the expected change in \(Y\) for a one-unit increase in \(X1\) controlling for \(X2\). See Semester 1 Week 5 Lab for a recap on interpretation.

• As presented in Table 3, the estimated value for the scientific knowledge slope was estimated to be approximately -0.08
• Results suggested that, holding age constant, for every one point increase in scientific knowledge quiz scores, attitude scores decreased by 0.08 points.
• Keeping in mind the valence of the question wording, this means that those who were more knowledgeable of science tended to be more favorable towards science – i.e. disagreeing with the statement.

Age

Hint - Interpretation

This represents the average marginal effect of \(X2\) on \(Y\), and can be interpreted as the expected change in \(Y\) for a one-unit increase in \(X2\) controlling for \(X1\). See Semester 1 Week 5 Lab for a recap on interpretation.

• As presented in Table 3, the estimated value for the age slope was estimated to be approximately 0.002
• Results suggested that, holding science knowledge constant, for every one year increase in age, attitude scores increased by 0.002 points.
• Keeping in mind the valence of the question wording, this means that older people tend to be less favorable towards science – i.e. agreeing with the statement.

Question 9

Obtain 95% confidence intervals for your bootstrapped model estimates.

Hint

The type = "perc" argument tells R to return the values that comprise 95% of all values in between them, i.e. the value with 2.5% of observations below it and the value with 2.5% of observations above it and 97.5% of observations below it.

Solution

You can use your preferred confidence level here, but by default this is 95%:

Confint(boot_mdl, level = 0.95, type = "perc")

Bootstrap percent confidence intervals

                      Estimate         2.5 %       97.5 %
(Intercept)        2.788248712  2.6877669516  2.888953606
science_knowledge -0.080276256 -0.0887889092 -0.071671501
age                0.002379446  0.0009551896  0.003710649

If you want to make it into a nice table, we can use kable():

Confint(boot_mdl, type = "perc") %>%
    kable(digits = 3, caption = 'Bootstrap 95% CIs') %>%
    kable_styling(full_width = FALSE)

Table 4: Bootstrap 95% CIs

	Estimate	2.5 %	97.5 %
(Intercept)	2.788	2.688	2.889
science_knowledge	-0.080	-0.089	-0.072
age	0.002	0.001	0.004

• We are 95% confident that the population intercept is between 2.68 and 2.90
  
• We are 95% confident that the population slope for science knowledge is between -0.09 and -0.07
  
• We are 95% confident that the population slope for age is between 0.001 and 0.004

Writing Up & Presenting Results

Question 10

Interpret the results from your bootstrapped model in the context of the research question.

Make reference to your key result table(s) and plot(s).

Hint

Make sure to include a decision in relation to your null hypothesis - based on the evidence, should you reject or fail to reject the null?

Solution

We used a subset of data from the 2005 Eurobarometer 63.1 survey to investigate whether there was an association between people’s scientific knowledge and their attitudes towards science and faith after accounting for their age.

To answer the research hypothesis, we fitted the following regression model: \[ \text{Attitude} = \beta_0 + \beta_1 \cdot \text{Science Knowledge} + \beta_2 \cdot \text{Age} + \epsilon \]

Which resulted in the following estimated regression coefficients for the original sample:

\[ \widehat{\text{Attitude}} = 2.8 -0.08 \cdot \text{Science Knowledge} + 0.0024 \cdot \text{Age} \]

The model did not satisfy the regression assumptions (see Figure 6) and for this reason we assessed statistical significance using the bootstrap approach with \(R = 1000\) resamples.

The 95% bootstrap confidence intervals are provided in Table 4. Results suggested that there was a negative and statistically significant association between attitudes towards science and faith and scientific knowledge after controlling for age \((\beta = -0.08, CI_{95}[-0.09, -0.07])\). Since the 95% CI did not contain zero, we rejected the null hypothesis as there was evidence of an association between peoples’ attitudes towards science and faith and their scientific knowledge after accounting for their age. Specifically, results suggested that for every additional quiz question people got correct, we expected their attitude score to be lower by about 0.08 points, holding age constant. In other words, respondents with greater scientific knowledge tended to be more favorable towards it, regardless of their age.

Compile Report

Compile Report

Knit your report to PDF, and check over your work. To do so, you should make sure:

• Only the output you want your reader to see is visible (e.g., do you want to hide your code?)
• Check that the tinytex package is installed
• Ensure that the ‘yaml’ (bit at the very top of your document) looks something like this:

---
title: "this is my report title"
author: "B1234506"
date: "07/09/2024"
output: bookdown::pdf_document2
---

What to do if you cannot knit to PDF

If you are having issues knitting directly to PDF, try the following:

• Knit to HTML file
  
• Open your HTML in a web-browser (e.g. Chrome, Firefox)
  
• Print to PDF (Ctrl+P, then choose to save to PDF)
  
• Open file to check formatting

Hiding Code and/or Output

Hiding R Code

To not show the code of an R code chunk, and only show the output, write:

```{r, echo=FALSE}
# code goes here
```

Hiding R Output

To show the code of an R code chunk, but hide the output, write:

```{r, results='hide'}
# code goes here
```

Hiding R Code and Output

To hide both code and output of an R code chunk, write:

```{r, include=FALSE}
# code goes here
```

Tinytex

You must make sure you have tinytex installed in R so that you can “Knit” your Rmd document to a PDF file:

install.packages("tinytex")
tinytex::install_tinytex()

Solution

You should end up with a PDF file. If you have followed the above instructions and still have issues with knitting, speak with a Tutor.