Loops, Conditionals and Variables: A Basic Simulation in R

May 3, 2011
By

(This article was first published on Psychwire » R, and kindly contributed to R-bloggers)

In this post, I’m going to go over some basics of using conditionals and loops in R. I’ll expand on the example I use here in future posts. The conditionals and loops will be used to create some dummy eye movement data.

Background

Before I get into the actual code itself, I should probably explain what eye movements are all about. It’s a pretty big topic but basically, the easy way to think about it is to consider that your eyes look at things in the environment that you’re interested in. Put in more scientific terminology, your eyes fixate (point to) objects or areas of the environment containing information that your brain and cognitive systems are trying to process in detail. This happens because, though you don’t realise it, the quality of the visual input from your eyes is actually very poor. You only have colour vision in the dead centre of your visual field (though never, ever notice it), and beyond the centre of your visual field, the input not only goes from being in colour to being in black-and-white, but the clarity and resolution drops off significantly as well.

The solution to the limited quality of visual input is to (1) utilise a load of systems that make you feel comfortable and safe, with everything neatly in colour and crystal clear and (2) to move your eyes around. A lot. You make 5-6 eye movements every second that you are awake, and, though you can of course have conscious control of them, most of the time, you let your eyes scoop up information in the outside world that is relevant so you can basically just get on with your life.

Your eyes are never truly still, though there are periods when they are still and information is taken in. These time periods are called fixations and are interesting because you take in information during fixations. The movements between fixations are called saccades (French for ‘jumps’), and, though you don’t realise it, you’re blind during these saccades. Saccades are short (around 60ms, though this depends on the task), and fixations are much longer (varies considerably, but here we’re talking about 200ms).

The Simulation

Here I’m going to simulate people looking at four different objects: a square, a circle, a star and a triangle. If you imagine a display is drawn out in front of a participant, and these four objects are present. The participant’s job is to locate a circle. Once they have done that, they press a button and the trial ends. However, if they don’t find the circle, they can also give up, but they won’t do that straight away.

Let’s begin!

The Code

We begin by creating a dataframe called fix_table. It has 10,000 rows, given by the seq function.

fix_table <-data.frame(seq(1:10000))

Next we set some default values and create some columns. Trial is the simulated trial number that we’re in. Object is the object the participant is looking at, be it the square, circle, star or triangle. Fix_index is the current index of the fixation – this gets reset to 1 at the start of each trial.

fix_table$trial = 0
fix_table$object = "null"
fix_table$fix_index = 0

Now we set some defaults before running the main loop of the code. Objects is the list of different objects presented in each trial to look at. Fix_index begins at 1 because of it being the first fixation in a trial. Trial_index starts at one for the first trial:

objects = list("circle", "square", "star", "triangle")
fix_index = 1
trial_index = 1

The Simulation Loop

Now for the actual loop that does the simulation itself. It’s a for loop that goes through each row of the fix_table dataframe, starting at 1 (the first row) and ending at the final row, determined by the number of rows function or nrow(fix_table).

for (row in 1:nrow(fix_table)) {
... code goes here...
}

So, what code do we want to go into the loop? We begin by setting the basic information for that row, updating the fixation index and trial index values, like this:

fix_table[row,"fix_index"] = fix_index
fix_index = fix_index + 1
fix_table[row, "trial"] = trial_index

Next we randomly sample one of the objects to be looked at by the participant:

current = sample(objects,1)[1]
fix_table[row, "object"] = current

Sample randomly selects 1 object from the objects list, and then gets assigned to current. We then update the dataframe called fix_table with the name of the current object being looked at.

After this, we need to decide whether a trial is going to end with the current fixation:

p_end <-rnorm(1, mean=1/fix_index, sd=0.3)
if (p_end>1 | current=="circle") { trial_end=TRUE }
else {trial_end=FALSE}

This is some made-up code that first of all creates a sort-of random probability value that the trial will end. The trial is more likely to end as more fixations are made, and we generate a normally distributed random number with mean of 1/fix_index and standard deviation of 0.3. If this value is greater than 1, the trial will end. Alternatively, as participants are searching for a circle, if they look at the circle, the trial will end. This is determined by the use of the or condtional, signified by the vertical pipe, |. Otherwise, the trial continues.

If the trial is set to end, we need to reset some important values for the new trial to begin. We do this via the following:

if (trial_end==TRUE){
trial_index = trial_index + 1
fix_index = 1
}

Simple!

The Full Code

Here we go:

fix_table <-data.frame(seq(1:100))

fix_table$trial = 0
fix_table$object = "null"
fix_table$fix_index = 0
objects = list("circle", "square", "star", "triangle")
fix_index = 1
trial_index = 1
for (row in 1:nrow(fix_table)) {
# add basics
fix_table[row,"fix_index"] = fix_index
fix_index = fix_index + 1
fix_table[row, "trial"] = trial_index
# decide which object we are on this time
current = sample(objects,1)[1]
fix_table[row, "object"] = current
# determine if trial ends!
p_end <-rnorm(1, mean=1/fix_index, sd=0.3)
if (p_end>1 | current=="circle") { trial_end=TRUE }
else {trial_end=FALSE}
# if the trial ends, reset values
if (trial_end==TRUE) {
trial_index = trial_index + 1
fix_index = 1   }
}

Finally, let’s make a histogram of how long it takes for a trial to end:

summary_table <-ddply(fix_table, c("trial"), 
summarise, max=max(fix_index))hist(summary_table$max)

Which gives the following:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

That’s it for now! More complex aspects will be added later.

To leave a comment for the author, please follow the link and comment on his blog: Psychwire » R.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...



If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Tags: , , , ,

Comments are closed.