General preparation

Read in code

source("rulexj_functions.R")

Set up exemplar matrix

columns correspond to cues, last column is criterion value
rows correspond to exemplars

exemplars <- matrix(
  ncol=5,byrow=TRUE,dimnames = list(NULL,c("c1","c2","c3","c4","crit")),
  data = c(
    0,0,0,1,23,
    0,0,1,0,25,
    0,1,0,0,30,
    0,1,0,1,43,
    1,0,0,0,35,
    1,0,1,1,70,
    1,1,0,1,68,
    1,1,1,0,63))

# display exemplars
exemplars

##      c1 c2 c3 c4 crit
## [1,]  0  0  0  1   23
## [2,]  0  0  1  0   25
## [3,]  0  1  0  0   30
## [4,]  0  1  0  1   43
## [5,]  1  0  0  0   35
## [6,]  1  0  1  1   70
## [7,]  1  1  0  1   68
## [8,]  1  1  1  0   63

Read in data

library(foreign)
dataset <- read.spss("exp1A.sav",to.data.frame = TRUE)

## re-encoding from UTF-8

# select only last three training blocks of 3 participants for demonstration
dataset <- subset(dataset, VpNr%in%1:3 & Phase == "training" & Block >=6)

# Code cues
dataset$c1 <- (grepl("strong",dataset$Pat1a))*1
dataset$c2 <- (grepl("strong",dataset$Pat1b))*1
dataset$c3 <- (grepl("strong",dataset$Pat1c))*1
dataset$c4 <- (grepl("strong",dataset$Pat1d))*1

Fit data

The function fit_rulexj can be used to fit the participant judgments. The following arguments should be provided:

dataset: dataset containing the to be fitted judgments and cue patterns
exemplars: matrix of the exemplars stored in memory
id_label: name of the column in the dataset containing the participant identifier
cue_labels: vector of the column names in the dataset containing the cue patterns
judgments_label: name of the column in the dataset containing the judgments of the participant
repetitions: number of times the fitting procedure should be repeated with new randomly generated start parameters. Note that when fitting the actual participant parameters for the article, repetitions was set to 50. The default value is 10 which was used for the recovery simulations.

In addition, a number of optional arguments can be specified for further customization:

report: boolean indicating whether the function should report on its progress, default is FALSE
alpha_bounds: lower and upper limit for the alpha parameter, defaults are c(.00001,.99999)
w_bounds: lower and upper limit for the weight parameters, defaults are c(0,100)
s_bounds: lower and upper limit for the attention parameter, defaults are c(.00,.50)
fit_rulexj_joint: boolean indicating whether the alpha parameter should be fit jointly with the other parameters or afterwards, default is TRUE

fitted_data <- fit_rulexj(
  dataset=dataset,
  exemplars=exemplars,
  id_label="VpNr",
  cue_labels=c("c1","c2","c3","c4"),
  judgments_label="Urteil",
  repetitions = 50,
  report=TRUE
)

## [1] "VpNr 3 finished"
## [1] "VpNr 1 finished"
## [1] "VpNr 2 finished"

fit_rulexj returns a data.frame with one row per participant. It contains the following variables:

First column contains the subject identifier
The results for the parameters when the rule and the exemplar module was fitted separately (with suffix _sep)
The results for the RulEx-J model where the rule and exemplar module was fitted together with the alpha parameter (with suffix _joint).
Note that this differs when fit_rulexj_joint is set to FALSE, as in this case only alpha_sep is fitted in addition to the module parameters.

round(fitted_data,2)

##   VpNr w0_sep w1_sep w2_sep w3_sep w4_sep SS_Rule_sep s_sep SS_Exemplar_sep
## 1    3  10.28  24.35  16.32  15.28  16.32      172.49  0.03           25.77
## 2    1  36.63   0.00   9.28  18.13   9.12     6010.13  0.49         8558.07
## 3    2  16.33  26.83  12.17  12.00  12.17      777.33  0.05          905.78
##   alpha_joint w0_joint w1_joint w2_joint w3_joint w4_joint s_joint
## 1        0.01     0.00     0.00    22.38    41.38    16.44    0.02
## 2        0.93    37.12     0.00     9.24    18.48     9.05    0.33
## 3        0.25    34.39    32.22     0.00     3.47     0.00    0.00
##   SS_RulExJ_joint
## 1           24.63
## 2         6004.58
## 3          639.97

Recovery simulation

Simulate data

The function simulate_rulexj can be used to randomly generate parameters for the RulEx-J model and simulate the resulting judgments for a given set of cues. The following arguments need to be provided:

n: the number of participants to simulate
cues: the cue patterns for which judgments should be derived
exemplars: the exemplars stored in memory

Besides, a number of optional arguments can be specified which already have reasonable defaults:

alpha_fix: an optional parameter indicating if alpha should be kept constant
alpha_bounds: a vector specifying the range of allowed alpha values
w_bounds: a vector specifying the range of allowed values for the weight parameters
s_bounds: a vector specifying the range of allowed attention parameter values
w_sum_max: ab optional numeric value - if specified, the sum of all w cannot be greater than this value

simulated_data <- simulate_rulexj(
  n=10,
  cues=exemplars[,1:4], # use all exemplar cue patterns as cues
  exemplars=exemplars,
  w_sum_max = 100
)

simulate_data returns a list of two data.frames.

parameters contains the data generating parameters:

round(simulated_data$parameters,2)

##       ID alpha    w0    w1    w2    w3    w4    s
##  [1,]  1  0.34 20.70  4.07 13.01 11.35 50.31 0.02
##  [2,]  2  0.73 20.05 13.85  6.54  4.33  9.01 0.33
##  [3,]  3  0.18 82.25  4.39  4.55  0.74  4.39 0.13
##  [4,]  4  0.06  9.77 10.86 17.01 23.61  9.15 0.33
##  [5,]  5  0.26  0.85  8.20 18.73 34.74 35.13 0.36
##  [6,]  6  0.48 12.30 14.59 23.44  3.57 35.56 0.16
##  [7,]  7  0.24 19.39 24.78  6.04 33.06  7.14 0.08
##  [8,]  8  0.66  9.31 37.96 12.32  1.59 13.93 0.37
##  [9,]  9  0.82  5.34 23.08  7.51  4.38 32.41 0.09
## [10,] 10  0.08 25.28  6.96 22.97  3.14  0.19 0.34

judgments contains the resulting judgments for the cue patterns:

head(simulated_data$judgments,10)

##       ID judgment c1 c2 c3 c4
##  [1,]  1 39.64086  0  0  0  1
##  [2,]  1 27.42509  0  0  1  0
##  [3,]  1 31.44722  0  1  0  0
##  [4,]  1 56.86884  0  1  0  1
##  [5,]  1 31.53405  1  0  0  0
##  [6,]  1 75.55815  1  0  1  1
##  [7,]  1 74.48464  1  1  0  1
##  [8,]  1 58.24953  1  1  1  0
##  [9,]  2 30.15348  0  0  0  1
## [10,]  2 26.46463  0  0  1  0

Fitting the simulated data

The function fit_rulexj can be used to fit the data in order to check whether the data generating parameters can be recovered.

recovered_data <- fit_rulexj(
  dataset=simulated_data$judgments,
  exemplars = exemplars
)

It returns a data.frame with the best fitting parameter values.

round(recovered_data,2)

##    ID w0_sep w1_sep w2_sep w3_sep w4_sep SS_Rule_sep s_sep SS_Exemplar_sep
## 1   1  13.89  17.68  15.15  13.72  27.90       20.88  0.07          575.88
## 2   2  21.31  14.17   7.14   5.51   9.02        1.19  0.38          299.40
## 3   3  26.17  18.95  12.45  11.32  12.69       27.84  0.20          497.10
## 4   4  24.23  14.63   9.12   9.44   8.94       14.37  0.33            2.50
## 5   5  19.73  12.47  10.81  14.95  15.25        7.36  0.27           86.60
## 6   6  14.21  17.83  18.16   8.33  24.22        9.88  0.12          478.26
## 7   7  13.81  23.84  13.06  18.94  13.50       25.83  0.09           97.17
## 8   8  15.28  29.67  10.79   3.71  11.93        1.47  0.13          252.76
## 9   9   6.59  23.11   8.84   6.12  29.41        1.36  0.08          899.40
## 10 10  25.29  14.11   9.63   8.10   8.16       13.35  0.35            4.63
##    alpha_joint w0_joint w1_joint w2_joint w3_joint w4_joint s_joint
## 1         0.34    20.70     4.07    13.01    11.35    50.31    0.02
## 2         0.73    20.05    13.85     6.54     4.33     9.01    0.33
## 3         0.18    82.25     4.39     4.55     0.74     4.39    0.13
## 4         0.06     9.77    10.86    17.01    23.61     9.15    0.33
## 5         0.26     0.85     8.20    18.73    34.74    35.13    0.36
## 6         0.48    12.30    14.59    23.44     3.57    35.56    0.16
## 7         0.24    19.39    24.78     6.04    33.06     7.14    0.08
## 8         0.66     9.31    37.96    12.32     1.59    13.93    0.37
## 9         0.82     5.34    23.08     7.51     4.38    32.41    0.09
## 10        0.08    25.28     6.96    22.97     3.14     0.19    0.34
##    SS_RulExJ_joint
## 1                0
## 2                0
## 3                0
## 4                0
## 5                0
## 6                0
## 7                0
## 8                0
## 9                0
## 10               0

Assessing the recovery performance

The recovery performance can be assessed with some basic functions in R.

# Plot of data generating and fitted alpha
plot(simulated_data$parameters[,"alpha"],recovered_data$alpha_joint,las=1)

# Summary of sums of squares of the RulEx-J model
summary(recovered_data$SS_RulExJ_joint)

##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## 3.200e-23 4.453e-20 1.247e-19 1.387e-19 1.652e-19 3.778e-19

# Summary of differences between data generating and fitted parameters:
summary(subset(simulated_data$parameters,select=c(alpha:s))-
          subset(recovered_data,select=c(alpha_joint:s_joint)),
        digits = 2)

##   alpha_joint          w0_joint           w1_joint           w2_joint       
##  Min.   :-5.8e-10   Min.   :-8.7e-08   Min.   :-1.6e-08   Min.   :-5.3e-09  
##  1st Qu.:-5.8e-11   1st Qu.:-5.0e-09   1st Qu.:-3.2e-09   1st Qu.:-2.8e-10  
##  Median : 1.5e-12   Median : 1.6e-10   Median :-1.7e-10   Median :-1.1e-10  
##  Mean   : 5.7e-12   Mean   :-7.9e-09   Mean   :-1.4e-09   Mean   : 4.3e-09  
##  3rd Qu.: 2.5e-11   3rd Qu.: 1.5e-09   3rd Qu.: 8.0e-10   3rd Qu.: 3.0e-09  
##  Max.   : 4.3e-10   Max.   : 1.5e-08   Max.   : 7.4e-09   Max.   : 3.8e-08  
##     w3_joint           w4_joint           s_joint        
##  Min.   :-1.5e-09   Min.   :-6.1e-09   Min.   :-4.7e-10  
##  1st Qu.:-6.0e-10   1st Qu.:-3.8e-10   1st Qu.:-3.4e-12  
##  Median :-2.5e-11   Median : 2.2e-10   Median : 4.5e-12  
##  Mean   : 7.3e-09   Mean   : 7.9e-09   Mean   :-4.2e-11  
##  3rd Qu.: 6.0e-10   3rd Qu.: 4.4e-09   3rd Qu.: 1.6e-11  
##  Max.   : 7.4e-08   Max.   : 7.5e-08   Max.   : 2.2e-10

RulEx-J demonstration