# General information: # 1) We will assume for the sake of simplicity that the data fulfill all requirements # for the use of t tests: Observations are randomly and independently drawn from # normally distributed populations of equal variance. # 2) The 95% CIs refer only to the individual mean (arelational CI) and cannot be used # directly to draw conclusions about differences between means (i.e., do not use the CIs # to compare different contrasts with each other). # Example of contrast analyses of one-way between-subject data # # Data similar to that in Howell (2013, table 13.2) who modelled it after a study by Eysenck (1974). # Number of recalled words in relation to five levels of encoding, defined according to # instructions to participants, and ordered from shallow to deep # (C = count letters, R = rhyme, A = think of appropriate adjective, # I = imagine visually, M = memorize for later recall). # # We will treat the data as a one-way between-subject design # with 20 participants per encoding condition # # Stefan Wiens & Mats E. Nilsson, 2016-06-23 # Contact: mats.nilsson@psychology.su.se # Clear everything ------------------------------------------------------------- rm(list = ls()) # clear memory graphics.off() # clear all plots cat("\014") # clear console (same as Ctrl-L in console) #------------------------------------------------------------------------------- # Set working directory (mydir) ------------------------------------------------ mydir <- 'c:/Users/MATNI/Dropbox/Statspaper with Mats/Supplement/R' setwd(mydir) # ------------------------------------------------------------------------------ # Load functions --------------------------------------------------------------- # The file 'functions_contrasts.r' (should be in mydir) contain code for these # functions: get_example_data(), get_contrast(), plot_contrast() source('functions_contrasts.r') # ------------------------------------------------------------------------------ # Get example data, taken from Howell, Table 13.2 ------------------------------ d <- get_example_data(design = 'between', way = 'one') # ------------------------------------------------------------------------------ # Screening of data ------------------------------------------------------------- summary(d) boxplot(d$recall ~ d$encode) aggregate(d$recall, list(d$encode), summary) # ------------------------------------------------------------------------------ # Define contrast-weight matrix------------------------------------------------- # Create matrix with five contrasts # If you run "contrast_result" below with "standardset = TRUE", # these contrasts will be converted to standard sets (i.e. absolute weights sum to 2) # If you run "contrast_result" below with "standardset = FALSE", # these contrasts will be computed with the given weights. coef <- matrix(c(-1, 1, 0, 0, 0, # 1. R - C -1, 0, 1, 0, 0, # 2. A - C 0, 0, -2, 1, 1, # 3. I+M - A -3, -3, 2, 2, 2, # 4. A+I+M - C + R -2, -1, 0, 1, 2), # 5. Linear trend nrow = 5, ncol = 5, byrow = TRUE) # Check if contrast sum to zero (within 1e-12 margin of error) if (abs(sum(apply(coef,1,sum))) > 1e-12) stop ('Not all contrasts sum to zero!') # Add row and column names to the contrast-weight matrix cnames <- c('R>C', 'A>C','I+M>A', 'A+I+M>C+R', 'Trend') rownames(coef) <- cnames colnames(coef) <- levels(d$encode) #------------------------------------------------------------------------------- # Calculate and plot contrast scores with CIs, by using get_contrast() ------------ # Below, the error = 'anova2' will pool the error only from groups that are # included in the contrast (i.e., from groups that have nonzero weights). # The error = 'anova' will pool the error from all groups # (even if they have zero weights). This is appropriate if all groups have # similar variances. # See below for an example. contrast_result <- apply(coef, 1, get_bcontrast, data = d$recall, group = d$encode, error = 'anova', standardset = FALSE) # If you run "contrast_result" with "standardset = TRUE", # these contrasts will be converted to standard sets (i.e. absolute weights sum to 2) # If you run "contrast_result" with "standardset = FALSE", # these contrasts will be computed with the given weights. # Transpose matrix using t(), if you prefer to have each contrast in a single row contrast_result <- t(contrast_result) print(contrast_result) plot_contrast(contrast_result) #------------------------------------------------------------------------------- # Example results: # # contrast_result <- apply(coef, 1, get_bcontrast, data = d$recall, # group = d$encode, error = 'anova', standardset = FALSE) # # Note that the dfs and MSerror are identical for all contrasts. The reason is that # with "error=anova", the error is pooled from all groups (including those with zero weights) # Note further that the contrast weights are not turned into standard sets (standardset = FALSE). # For example, if you entered [+2 -2 0 0 0], this captures the difference between the first # and second mean, but this difference is multiplied by 2. # # Score CI_LL CI_UL Tobs df p SE MSerror C R A I M # R>C 0.40 -1.542259 2.342259 0.408854 95 6.835671e-01 0.9783445 9.571579 -1 1 0 0 0 # A>C 6.05 4.107741 7.992259 6.183916 95 1.563474e-08 0.9783445 9.571579 -1 0 1 0 0 # I+M>A 5.40 2.035908 8.764092 3.186701 95 1.947994e-03 1.6945423 9.571579 0 0 -2 1 1 # A+I+M>C+R 45.90 38.377662 53.422338 12.113657 95 5.626533e-21 3.7891118 9.571579 -3 -3 2 2 2 # Trend 25.95 21.606976 30.293024 11.862072 95 1.883710e-20 2.1876447 9.571579 -2 -1 0 1 2