Week 8

Week 8 Lab Activity / Homework

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.2     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── 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
## Loading required package: Matrix
## 
## 
## Attaching package: 'Matrix'
## 
## 
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
## 
## 
## 
## Attaching package: 'psych'
## 
## 
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
## 
## 
## 
## Attaching package: 'MASS'
## 
## 
## The following object is masked from 'package:dplyr':
## 
##     select
## 
## 
## 
## arm (Version 1.13-1, built: 2022-8-25)
## 
## 
## Working directory is C:/Users/brand/Desktop/UCR/PhD CCN/PSYC213_Spring2023/Lab/Week 8
## 
## 
## 
## Attaching package: 'arm'
## 
## 
## The following objects are masked from 'package:psych':
## 
##     logit, rescale, sim

# Lets start by reading in the Radon data  

load("radon.Rdata")

# Now Read in the data and create factor variables for Haptic Search

hs_data <- read.csv("MLM_haptic_example_data.csv", header = TRUE)
hs_data <- as_tibble(hs_data) %>% mutate(Sub_Num = factor(Subject_Num))    
hs_data <- hs_data %>% mutate(Ratio = factor(Relative_Size))
hs_data <- hs_data %>% mutate(Num_Dis = factor(Distractor_Number))

We need to make the appropraite changes to the radon data

# start by creating a county sample size column 
county_n <- tibble(table(radon$county))

county_n <- county_n %>% mutate(county = 1:85) 
county_n <- county_n %>% mutate(N = table(radon$county)) 

#put all data frames into list
df_list <- list(radon, radon_county, county_n)

#merge all data frames in list
radon_full <-df_list %>% reduce(full_join, by='county')

Varying Intercept Model with 1 predictor (Chapter 12)

Just as a refresher lets go over the main model from chapter 12

Fit the varying intercept with 1 predictor model using lmer(). Fit to both data sets

M1<-lmer(y ~ 1 + x +(1 | county), data = radon_full)
M1

## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ 1 + x + (1 | county)
##    Data: radon_full
## REML criterion at convergence: 2171.305
## Random effects:
##  Groups   Name        Std.Dev.
##  county   (Intercept) 0.3282  
##  Residual             0.7556  
## Number of obs: 919, groups:  county, 85
## Fixed Effects:
## (Intercept)            x  
##       1.462       -0.693

list(M1)

## [[1]]
## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ 1 + x + (1 | county)
##    Data: radon_full
## REML criterion at convergence: 2171.305
## Random effects:
##  Groups   Name        Std.Dev.
##  county   (Intercept) 0.3282  
##  Residual             0.7556  
## Number of obs: 919, groups:  county, 85
## Fixed Effects:
## (Intercept)            x  
##       1.462       -0.693

hs.m.0<-lmer(Times ~ 1 + Relative_Size + (1 | Distractor_Lengths), data = hs_data)
display(hs.m.0)

## lmer(formula = Times ~ 1 + Relative_Size + (1 | Distractor_Lengths), 
##     data = hs_data)
##               coef.est coef.se
## (Intercept)   17.58     1.23  
## Relative_Size -2.93     0.91  
## 
## Error terms:
##  Groups             Name        Std.Dev.
##  Distractor_Lengths (Intercept) 2.06    
##  Residual                       9.74    
## ---
## number of obs: 3370, groups: Distractor_Lengths, 6
## AIC = 24930.8, DIC = 24928.3
## deviance = 24925.6

#icc1<- 0.3282^2 / (0.3282^2 + 0.7556^2)
#icc1

Question 1: Describe the results for each model given the data. What is the effect of group level information on individual level data?

Answer: When holding county radon level constant, the average house radon decreases by 0.69 from the basement to the first floor.

Question 2: Pick one of the models and using the `ranef()` and the `fixef()` functions provide the coefficient estimates for at least 3 of the groups.

Answer:

fixef(hs.m.0)

##   (Intercept) Relative_Size 
##     17.584509     -2.927372

ranef(hs.m.0)

## $Distractor_Lengths
##      (Intercept)
## 0.25   0.7744182
## 0.5   -2.1478971
## 0.75   2.9210016
## 1.25   0.5752139
## 1.5   -1.1968970
## 1.75  -0.9258397
## 
## with conditional variances for "Distractor_Lengths"

fixef(hs.m.0)[[1]] + ranef(hs.m.0)[[1]]

##      (Intercept)
## 0.25    18.35893
## 0.5     15.43661
## 0.75    20.50551
## 1.25    18.15972
## 1.5     16.38761
## 1.75    16.65867

coef(hs.m.0)

## $Distractor_Lengths
##      (Intercept) Relative_Size
## 0.25    18.35893     -2.927372
## 0.5     15.43661     -2.927372
## 0.75    20.50551     -2.927372
## 1.25    18.15972     -2.927372
## 1.5     16.38761     -2.927372
## 1.75    16.65867     -2.927372
## 
## attr(,"class")
## [1] "coef.mer"

County 1 -0.27009750 + 1.4615979 County 2 -0.53395107 + 1.4615979 County 3 0.01761646 + 1.4615979

Chapter 13 Varying intercepts and varying slopes

Question 3: Fit a model with varying intercepts and varying slopes that does not include a group level predictor.

M2<-lmer(y~x+(1+x|county),data=radon_full)
display(M2)

## lmer(formula = y ~ x + (1 + x | county), data = radon_full)
##             coef.est coef.se
## (Intercept)  1.46     0.05  
## x           -0.68     0.09  
## 
## Error terms:
##  Groups   Name        Std.Dev. Corr  
##  county   (Intercept) 0.35           
##           x           0.34     -0.34 
##  Residual             0.75           
## ---
## number of obs: 919, groups: county, 85
## AIC = 2180.3, DIC = 2153.9
## deviance = 2161.1

(a) Interpret your results

Answer: The floor radon measurement, x, looks to be significant. The expected change in the radon measurement for a one-unit increase in x, holding other variables constant, is -.68.the expected value of the response variable (y) when the predictor variable (x) is 0, is 1.46.

The correlation between the random intercepts and slopes is -0.34, indicating a negative correlation between the intercepts and slopes across counties. The residual standard deviation, representing the unexplained variability within each county after accounting for the fixed and random effects, is 0.75.

hs.m.1<-lmer(Times ~ Relative_Size + (1 + Relative_Size| Experiment_Number), data = hs_data)

## boundary (singular) fit: see help('isSingular')

display(hs.m.1)

## lmer(formula = Times ~ Relative_Size + (1 + Relative_Size | Experiment_Number), 
##     data = hs_data)
##               coef.est coef.se
## (Intercept)   17.54     2.68  
## Relative_Size -2.84     2.06  
## 
## Error terms:
##  Groups            Name          Std.Dev. Corr  
##  Experiment_Number (Intercept)   3.78           
##                    Relative_Size 2.91     -1.00 
##  Residual                        9.67           
## ---
## number of obs: 3370, groups: Experiment_Number, 2
## AIC = 24878.6, DIC = 24869.2
## deviance = 24867.9

summary(hs.m.1)

## Linear mixed model fit by REML ['lmerMod']
## Formula: Times ~ Relative_Size + (1 + Relative_Size | Experiment_Number)
##    Data: hs_data
## 
## REML criterion at convergence: 24866.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.8027 -0.7220 -0.1990  0.4189  4.6676 
## 
## Random effects:
##  Groups            Name          Variance Std.Dev. Corr 
##  Experiment_Number (Intercept)   14.252   3.775         
##                    Relative_Size  8.465   2.910    -1.00
##  Residual                        93.602   9.675         
## Number of obs: 3370, groups:  Experiment_Number, 2
## 
## Fixed effects:
##               Estimate Std. Error t value
## (Intercept)     17.543      2.680   6.545
## Relative_Size   -2.839      2.065  -1.375
## 
## Correlation of Fixed Effects:
##             (Intr)
## Relative_Sz -0.998
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')

Got “boundary (singular) fit: see help(‘isSingular’)” error when trying to fit the last model

(b) Interpret your results

Answer: As the relative size difference between the target and distractors increases, the reaction time to find the target decreases. There are really small standard error coefficients compared to the relative size and experiment number coefficients. The residual standard deviation is 9.67.

Question 4: Fit a model with varying intercepts and varying slopes that includes a group level predictor.

M3<-lmer(y~x+u+ x:u +(1+x|county),data=radon_full,start=list(theta=c(1,2,3)))
M3

## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ x + u + x:u + (1 + x | county)
##    Data: radon_full
## REML criterion at convergence: 2126.579
## Random effects:
##  Groups   Name        Std.Dev. Corr
##  county   (Intercept) 0.1244       
##           x           0.3070   0.41
##  Residual             0.7495       
## Number of obs: 919, groups:  county, 85
## Fixed Effects:
## (Intercept)            x            u          x:u  
##      1.4686      -0.6709       0.8081      -0.4195

display((M3))

## lmer(formula = y ~ x + u + x:u + (1 + x | county), data = radon_full, 
##     start = list(theta = c(1, 2, 3)))
##             coef.est coef.se
## (Intercept)  1.47     0.04  
## x           -0.67     0.08  
## u            0.81     0.09  
## x:u         -0.42     0.23  
## 
## Error terms:
##  Groups   Name        Std.Dev. Corr 
##  county   (Intercept) 0.12          
##           x           0.31     0.41 
##  Residual             0.75          
## ---
## number of obs: 919, groups: county, 85
## AIC = 2142.6, DIC = 2101.9
## deviance = 2114.2

hs.m.2<-lmer(Times ~ Relative_Size + Distractor_Lengths + (1 + Relative_Size | Experiment_Number), data=hs_data)

## boundary (singular) fit: see help('isSingular')

display(hs.m.2)

## lmer(formula = Times ~ Relative_Size + Distractor_Lengths + (1 + 
##     Relative_Size | Experiment_Number), data = hs_data)
##                    coef.est coef.se
## (Intercept)        19.41     2.73  
## Relative_Size      -3.34     2.07  
## Distractor_Lengths -1.37     0.40  
## 
## Error terms:
##  Groups            Name          Std.Dev. Corr  
##  Experiment_Number (Intercept)   3.77           
##                    Relative_Size 2.91     -1.00 
##  Residual                        9.66           
## ---
## number of obs: 3370, groups: Experiment_Number, 2
## AIC = 24868.9, DIC = 24857.5
## deviance = 24856.2

Got “boundary (singular) fit: see help(‘isSingular’)” error when trying to fit the last model

Week 8

Brandon A. Nolasco

2023-06-16

Week 8 Lab Activity / Homework

Varying Intercept Model with 1 predictor (Chapter 12)

Fit the varying intercept with 1 predictor model using lmer(). Fit to both data sets

Question 1: Describe the results for each model given the data. What is the effect of group level information on individual level data?

Question 2: Pick one of the models and using the `ranef()` and the `fixef()` functions provide the coefficient estimates for at least 3 of the groups.

Chapter 13 Varying intercepts and varying slopes

Question 3: Fit a model with varying intercepts and varying slopes that does not include a group level predictor.

(a) Interpret your results

(b) Interpret your results

Question 4: Fit a model with varying intercepts and varying slopes that includes a group level predictor.

Week 8

Brandon A. Nolasco

2023-06-16

Week 8 Lab Activity / Homework

Varying Intercept Model with 1 predictor (Chapter 12)

Fit the varying intercept with 1 predictor model using lmer(). Fit to both data sets

Question 1: Describe the results for each model given the data. What is the effect of group level information on individual level data?

Question 2: Pick one of the models and using the ranef() and the fixef() functions provide the coefficient estimates for at least 3 of the groups.

Chapter 13 Varying intercepts and varying slopes

Question 3: Fit a model with varying intercepts and varying slopes that does not include a group level predictor.

(a) Interpret your results

(b) Interpret your results

Question 4: Fit a model with varying intercepts and varying slopes that includes a group level predictor.

Question 2: Pick one of the models and using the `ranef()` and the `fixef()` functions provide the coefficient estimates for at least 3 of the groups.