In this practical you’ll practice “Programming in R” with the base and purrr
-package functions. Although we sort of have been programming already, we always used the tidyverse
way to do it. While there is nothing against doing it this way, in fact we chose to teach R with the tidyverse
because we think that it’s an easy and very elegant way of doing it, when you work with R in the future, you will almost certainly find yourself in situations were you have to understand the base R way of programming. This may, for example, be the case when you get code from someone, or work with someone who is not familiar with the tidyverse
. Furthermore, the concepts you learn in this session are the same in most programming languages (of course the exact implementation will vary over languages). So if R is your first programming language, learning about and understanding these concepts will greatly facilitate the understanding and speed of learning of other programming languages. Last but no least, the concepts we touch upon here taken together will allow you to code more efficiently and elegantly, which may save you quite some time. But enough said; let’s get started!
By the end of this practical you will know how to:
if
, else if
, and else
)for
- and while
loops, apply
family, map
family)Package | Installation |
---|---|
tidyverse |
install.packages("tidyverse") |
yarrr |
install.packages("yarrr") |
sendmailR |
install.packages("sendmailR") |
skimr |
install.packages("skimr") |
Conditional Statements
Function | Description |
---|---|
if (){} else if (){} else {} |
Conditional execution of code blocks |
ifelse() |
Vectorized conditional element selection |
switch() |
Select one of a list of alternatives |
Functions
Function | Description |
---|---|
function(){} |
Define a new function |
stop(), warning(), message() |
Print messages to the user, in the case of stop() break execution of the function |
stopifnot() |
Ensure the Truth of R Expressions. Break if any is FALSE |
Iteration
Function | Description |
---|---|
for (){} |
Repeatedly execute block of code for each element of a sequence |
apply(), sapply(), lapply() |
Repeatedly execute block of code for each element of a sequence |
map(), map_dbl(), map_int() |
Repeatedly execute block of code for each element of a sequence |
map2(), map2_dbl(), map2_int() |
Repeatedly execute block of code for each element of two sequences |
pmap(), pmap_dbl(), pmap_int() |
Repeatedly execute block of code for each element of multiple (i.e. a list of) sequences |
## Conditional Statements ------------
# Conditional statements are used to execute a block of code only if a given
# expression evaluates to TRUE
num <- 5
if (num < 8){
print("This expression evaluated to TRUE, thus num is smaller than 8.")
}
## Write your own function -----------
# First arguments are specified
weighted_mean <- function(vls, # has no default. This means it MUST be
# specified in a function call.
wghts = 1 / length(vls)){ # Has a default, i.e.
# it could be left un-
# specified and would then
# take the default.
if (!is.numeric(vls) || !is.numeric(wghts)){ # Check whether conditions are met
# and throw an error if they aren't
stop("Arguments must be numeric.")
}
if (sum(wghts) == 0){
return(0) # if the sum of the weights is 0, in this case we want it to return
# 0. If this condition is TRUE, all following code is NOT executed,
# that is the function will be exited early.
}
# Code block to be executed if the function is called. The result of the
# last line is returned
sum(values * weights) / sum(weights)
}
## Iteration --------------------
## for loops
# example 1:
# the following three loops do exactly the same. You may encounter all three
# ways of specifying the loops, but the latter two are more likely.
for (i in c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10)){
print(i) # note that within loops print statements have to be used explicitly
}
for (i in 1:10){
print(i) # note that within loops print statements have to be used explicitly
}
for (i in seq_len(10)){
print(i) # note that within loops print statements have to be used explicitly
}
# By the way, here's the difference between seq_len, and seq_along:
seq_len(5) # creates a vector from 1 to the number specified
seq_along(letters[1:5]) # creates a numeric sequence of the positions of the specified
# vector
# you could also use 1:5, instead of seq_len(5), but note that while
x <- vector("double", length = 0)
1:length(x)
# yields probably not the behavior you'd want in a loop, while
seq_along(x)
# will lead to the desired behavior
# example 2:
x <- vector("double", length = 10) # Output part
for (i in seq_along(x)){ # Sequence part
x[i] <- i # Body part
}
x
# example 3:
## while loops
i <- 1
while(i < 5){
print(i)
i <- i + 1
}
## apply
# create a vector of variable names from the variables you want to call the same
# function on
vars <- c("age", "height", "weight", "tattoos", "parrots", "grogg")
# compute the mean of these variables. Note that the second argument specifies,
# whether the mean should be taken over columns (2), or rows (1)
apply(pirates[, vars], 2, mean)
# you can also use a simpler apply version called sapply, where it automatically
# computes the mean column wise
sapply(pirates[, vars], mean)
## map
# the map equivalent is
map_dbl(pirates[, vars], mean)
# remember that in the map functions the suffix (in this case _dbl) stands for the
# type of vector that should be returned. If this doesn't match, e.g. if a character
# vector is created in map_dbl, it will throw an error. map without suffix will
# return a list.
File | Rows | Columns | Description | Source |
---|---|---|---|---|
pirates.csv |
1000 | 17 | A collection of information on 1000 pirates. Variables include the pirates’ sex, age, number of tattoos and number of parrots. | The dataset is from the yarrr package. You should already have it in your 1_Data folder, if not you can download it here. |
Open your R project. It should already have the folders 0_Data
and 1_Code
. Make sure that all of the data files listed above are contained in the folder
Open a new R script and save it as a new file called Programming_practical.R
in the 2_Code
folder. At the top of the script, using comments, write your name and the date. The, load the set of packages listed above with library()
.
For this practical, we’ll use the pirates.csv
data. This dataset contains information such as sex, age, and number of tattoos of 1000 pirates. Using the following template, load the data into R and store it as a new object called pirates
.
pirates <- read_csv(file = "XX/XXX")
print()
, summary()
, head()
and skim()
, explore the data to make sure it was loaded correctly.You have already used conditional statements, mostly in the form of the case_when()
function. Now let’s look at other functions. Note that for now the example may seem a little arbitrary. That is because usually conditional statements are used in loops or functions. We’ll bring everything together towards the end of this practical.
sample
function to draw one sample from {H, T} (i.e. a set, or in R a vector, containing “H” for head, and “T” for tail). If “H” comes up, print “Head”, if “T” comes up, print “Tail”.coin_toss <- sample(c("XX", "XX"), size = XX)
if (coin_toss == "XX"){
print("XX")
} else {
print("XX")
}
sample
function a bit. Check the help menu if you don’t know how to do this). Do this, and if ALL tosses came up heads, print “Wow, that’s ten heads in a row!”, if there is at least one head in there, print “The coin showed X times heads, and Y times tails”, where X and Y are the respective numbers of heads and tails (use the paste
function to print and the sum
function to get the correct numbers, check the help menus if you don’t know how to do this), and print “Wow, that’s ten tails in a row!”, if ALL tosses came up tails.coin_tosses <- sample(c("XX", "XX"),
size = XX,
replace = TRUE)
if (all(coin_tosses == "XX")){
print("Wow, that's ten heads in a row!")
} else if (any(coin_tosses == "XX")) { # This works because if-else statements
# are executed in the order they are written,
# so because of the first check it is safe
# to assume, that not all are heads.
paste("The coin showed",
sum(coin_tosses == "XX"), "times heads, and",
sum(coin_tosses == "XX"), "times tails.")
} else {
print("XXX")
}
Most of the time, conditional statements are either used in combination with functions or loops. We’ll start by looking at functions. Remember that a function has two parts, one in which you specify the arguments, and one in which you then write the code that is executed. The last line that is executed determines what the function returns. You can make this explicit with a return()
statement, but it is good practice to only use this if you want to exit the function early, i.e. before all code is executed. If you don’t remember how to specify functions, check again above in the examples.
CheckGuestList()
that takes a guest’s name Guest
, checks their name against a list of invitations Invitations
and returns one of two messages depending on whether or not the guest was on the list!Guest
is on the list, print “Welcome Guest! We have been waiting for you!”Guest
is NOT on the list, print “Who are you?! You can’t come in!!”CheckGuestList <- function(Guest,
Invitations = c("Donald", "Angela", "Justin")) {
if(Guest %in% Invitations) {
paste("XXX", Guest, "XX")
} else {
paste("XX")
}
}
Try your function with the arguments Guest = "Donald"
, Invitations = c("Angela", "Justin")
CheckColors()
that takes an argument color
and returns a sentence telling you whether or not color
is in fact a color in R.
CheckColor <- function(XX) {
if(XX %in% colors()) {
paste("Yes!", XX, "is a color in R!")
} else {
paste("Nope. Try again")
}
}
higher_lower_game
that takes a number guess
and compares it a number thinking
…higher_lower_game <- function(XX, XX) {
if(XX < XX) {return(paste("My number is Higher!"))}
if(XX > XX) {return(XX))}
if(XX == XX) {XX}
}
Test your function by first storing a random number between 1 and 100 as thinking
. Then test your function! How many guesses do you need to get it right?
Write your own function named check_logical
that tests whether a given logical object is TRUE
or FALSE
, and returns “Argument is TRUE”, and “Argument is FALSE”, respectively.
Use your check_logical
function to see whether ALL pirates are larger than 1.6 meters.
One important principle when programming is the Don’t Repeat Yourself (DRY) principle. That is, if you have written a block of code several times, with only slight changes in it, it is time to write a function for this. For example:
# Compute the absolute difference between means of a variable
# separated by another variable
# for tattoos
abs(mean(pirates$tattoos[pirates$sex == "male"]) -
mean(pirates$tattoos[pirates$sex == "female"]))
# for age
abs(mean(pirates$age[pirates$sex == "male"]) -
mean(pirates$age[pirates$sex == "female"]))
# for height
abs(mean(pirates$height[pirates$sex == "male"]) -
mean(pirates$height[pirates$sex == "female"]))
Each time you copy the code and change a small part, there is a small chance of making a mistake, for example by forgetting to change the name of one variable. Also imagine having coded the above example for 10 variables and then you realize that you should have taken the male
and other
group together to see the mean difference between this joint group and females
. Then you would have to go through each line again to manipulate it.
Now write a function called abs_mean_diff
, that computes the absolute mean difference between two groups. It should take a vector mean_vec
with data to compute the mean on, a vector cat_vec
with data to separate the first, and a list sep_list
with two vectors of categories to separate for (note that when you want to check whether a number of values also occurs in another vector, you may want to use %in%
instead of ==
for this logical expression or check). Then use this function to compute the same mean differences as in the code above.
abs_mean_diff
function from task C8 by adding a fourth argument met
(short for method), where you pass the function name as an argument WITHOUT quotation marks. Set the default value to this argument to mean
. Name the function abs_diff
. Then compute the mean and the median absolute difference of tattoos for male and female pirates.Passing functions as arguments can be very handy, especially when you can combine this with apply
or map
function so you can enter a list of functions.
...
. Now add the ...
as last argument in the abs_diff
function definition from the previous task, and also add it to the met()
function inside the abs_diff
function. Now you can call functions with abs_diff
that take more than just the numeric vector as argument, for example you can compute a trimmed mean with mean(x, trim = .2)
for the 20% trimmed mean. Compute the differences of the 20% trimmed mean with the modified abs_diff()
function (i.e. the outer 20% of data on each side are discarded, which avoids outliers).Iteration is an extremely important feature in programming. It let’s you execute a code block for a large number of objects, with only having to write the code block once. The iteration methods we will look at in the next few tasks are for
and while
loops. Later we will also look at apply
and map
functions.
Remember that a for
loop executes a block of code as many times as you specify, whereas a while
loop does the same thing but does not stop until a prespecified condition is met. If you don’t remember how to specify a loop, revisit the examples in the beginning of this practical.
for (val in vector){}
, instead of using either 1: length(vector)
, or seq_along(vector)
(did you see that we used val
instead of i
when defining the sequence in the loop. You can enter any name you want there, just make sure to refer to this name in the body of the loop). Loop over the following vector containing variable names of the pirates dataset c("age", "height", "tattoos", "beard.length")
, and compute the mean of each variable and store the value in a named vector (i.e. each value should be named to identify it, to do this use names
) called pirate_means
.vars <- c("XX", "XX", "XX", "XX")
pirate_means <- vector("double", length(XXX))
names(pirate_means) <- vars
for (nms in vars){
pirate_means[XXX] <- mean(pirates[[XXX]])
}
pirate_means
Use a loop to rewrite the code where you called abs_diff
multiple times in task C8. Remember that in a loop you have to use print
explicitly.
You want to check whether the college a pirate went to, is related to the number of groggs a pirate drinks a day on average. To get rid of the potentially influential factors of age, sex, and weight, you want to match the two subsamples of 100 pirates each, such that you have a maximum difference in the number of females between the two groups of 12, a maximum mean age difference of 10 years, a maximum difference in the standard deviation of age of 8 years, a maximum mean weight difference of 8 kilograms, and a maximum difference in the standard deviation of weight of 7 kilograms. Use a while
loop where you take a subsample of 100 pirates from each college, and use the above parameters of this sample for the check. Note: don’t bother about doing the analysis afterwards, the exercise is just about the looping part. Also the difference values for the check are arbitrary. If it takes too long for the loop to run just use larger values, for the exercise the precision doesn’t have to be high.
# prepare variable such that at least one iteration takes place
sex_diff <- Inf
age_diff_mean <- Inf
age_diff_sd <- Inf
weight_diff_mean <- Inf
weight_diff_sd <- Inf
while (XX > XX || XX > XX || v > XX ||
XX > XX || XX > XX){
CCCC_pirates <- pirates %>%
filter(college == "XXX") %>%
sample_n(XXX, replace = FALSE)
JSSFP_pirates <- pirates %>%
filter(college == "XXX") %>%
sample_n(XXX, replace = FALSE)
sex_diff <- abs(sum(CCCC_pirates$sex == "XXX") -
sum(JSSFP_pirates$sex == "XXX"))
age_diff_mean <- abs(mean(XXX) - mean(XXX))
age_diff_sd <- abs(sd(XXX) - sd(XXX))
weight_diff_mean <- abs(mean(XXX) - mean(XXX))
weight_diff_sd <- abs(sd(XXX) - sd(XXX))
}
# Here you would proceed with the analysis
apply
and map
functions
The apply
and map
function families are quite similar. The map
functions from the purrr
package were introduced as a more consistent and less ambiguous alternative to the apply
functions. Both function families let you iterate over objects and call functions on the subsets in each iteration. So basically you can use these functions instead of loops. In general, the basic structure is that you first specify the object over which you want to iterate, then the function to call in each iteration and then additional arguments that function might take. Depending on which of the functions you use, there will be additional arguments to specify.
for
loop you’ve used in task D2, to run the abs_diff
function several times? Rewrite this by replacing the for
loop once with sapply
and once with a map
function (figure out which map
function would fit best here).sapply(pirates[[c("XX", "XX", "XX")]], XXX, cat_vec = pirates$sex,
sep_list = list("male", "female"))
map_dbl(pirates[[c("XX", "XX", "XX")]], XXX, cat_vec =pirates$sex,
sep_list = list("male", "female"))
Challenge. Now there would be a more elegant way of programming the abs_diff
function such that you don’t have to use loops or apply
/ map
functions outside of it. Redefine the abs_diff
function and use map_dbl
, such that you can enter multiple columns to the mean_vec
argument, instead of only one. For this you can use an anonymous function in map_dbl
. This means that the function is only defined inside the map_dbl
call and does not exist outside of it. The syntax for a self programmed anonymous mean
function would be map_dbl(data, function(x){sum(x) / length(x)})
. That is instead of passing a function name to map_dbl
, you pass a function definition function(argument){do something}
. Remember that you can pass functions as arguments.
There is another map
function called map2
(again with the suffixes for the different types of vectors returned), that let’s you iterate over two objects simultaneously. Use the appropriate map2
function and to draw 10 samples from three normal distributions with means 5, 9, 12, and sds, 4, 2, 7, and save them as an object (note: the function to draw samples from a normal distribution is rnorm(n, mean, sd)
).
Now use map
to plot each of the samples drawn. Note: use par(mfrow = c(1, 3))
to draw three plots in a line in the same canvas (i.e., the 1 determines the number of rows of plots and the 3 determines the number of columns of plots. This way you can create multiple-plot panels).
par(mfrow = c(1, 3))
map(XX, plot, type = "b", ylab = "XX") # uses the base plot() function
As you can see, besides the plots you’ve got NULL
printed in your console each time you called plot
. This is because map returns whatever the function called returns and in case of plot NULL
is created as output (besides the actual plot that is drawn). To take care of this problem, there is a special family of purrr
functions called walk
, that are used for functions like plot that create something on the fly and only have a console output as undesired side product. Create the plots again using walk()
instead of map()
.
Redo task D6 but this time don’t draw 10 samples per distribution but draw 10, 15, and 20 samples, respectively. To iterate over an unlimited number of objects simultaneously, use the pmap
function. Check the help menu if you’re not sure how this works.
The next two tasks contain elements of all parts we’ve covered in this lesson. You will write a function to do a short simulation. Briefly some background: In statistics, we often accept a hypothesis if the p-value, i.e. the probability of the found-, or more extreme data under the H0 hypothesis that assumes that there is no difference between two distributions, is smaller or equal to 5%. This is tested with statistical tests, of which we will use the t-test (r function t.test(group1, group2)
). We will draw two random samples from the same distribution and run a t.test to check whether they are significantly different. If they are not, we will resample from the distribution and run the check once more. We will then see how strong the inflation of the alpha error is, with one resample of size x.
We will also encounter a new cool feature: you can send mails from R. For example, if you are running a simulation that takes a lot of time and you don’t want to check again and again whether it’s completed. So you can tell R to send you a mail when it’s done, and if you like, even attach the data or some core output. The package to do this is sendmailR
, with the function sendmail
.
sim_false_positive
, that takes the following arguments: n (the number of samples initially drawn per distribution), m (the mean of the distribution from which samples are drawn), sd (the standard deviation of the distribution from which samples are drawn), resample (the sample size of the resample), iter (number of iterations, i.e. number of times the whole procedure is repeated. The higher the number, the more precise the estimate), mail (a mail address from which to which you want to send your result), mail_subject (the subject of the mail). The start of your function should look like this (we’ve also included the code for the mail function to make it a little easier):## sim_false_positive functoin
sim_false_positive <- function(n = 30,
m = 5,
sd = 2,
resample = 10,
iter = 1000,
mail = NULL,
mail_subject = "Your results"){
# Here goes your code for the simulation. It should return an object called
# false_positive, containing the proportion of significant results (i.e. with
# a p-value smaller or equal to .05)
# Here is the part that sends the mail
if (!is.null(mail)){
# mail addresses and subject
from <- mail
to <- mail
subject <- mail_subject
# The body of the mail. Costumize it to your liking
body <- list(paste("Here are your results:", paste0("Resample was ",
resample, "."),
paste0("The initial samplesize was ", n,
", and you've drawn ", iter, " samples."),
paste0("The proportion of false positives was",
false_positive, ".")))
# Function that acutally sends the mail. Note that not all mail addresses
# will work.
sendmail(from, to, subject, body,
control=list(smtpServer="ASPMX.L.GOOGLE.COM"))
}
# return proportion of false positives
false_positive
}
Now program the middle part to do the sampling and iteration.
Now that you’ve created your function, use one of the map
functions to execute the function several times with different sample sizes and resample sizes. Note that if you only change resample, you will have to create an anonymous function, because otherwise the resample size will be passed to the first argument of sim_false_positive
, which is n
.
Program a function called boxes_game()
that let’s you play a decision making game. Here is how the game works: There is a room containing 100 boxes from which one contains a bomb. You can open as many boxes as you like, each opened box will earn you 1 point. However if you open the box containing the bomb, you will earn no points. Here’s how you can play the boxes game in R. First, create the room as an object room_100 which contains a vector with 99 values of 1 (representing 1 point) and one value of negative infinity (-Inf) which represents the bomb.
# This vector represents the room of 100 boxes
room_100 <- c(rep(10, 99), -Inf)
Now program the function with the arguments open
for the number of boxes to open, and room
to which you pass the room_100
vector created above. If no boxes are opened, print “You didn’t open any boxes! You earned nothing but are still alive”. If boxes were opened and one of them contained the bomb, print “You’re dead!!! You opened XXX boxes and got the bomb!!!”. If none of the boxes contained the bomb, print “Congratulations!!! You opened XXX boxes and earned XXX Points! Don’t you want to play again? :)”.
The readline()
function let’s you ask the user for an input and save it as an object. Redefine the boxes_game()
function such that you no longer need the open
argument but ask for a userinput instead.
Extend the boxes_game()
function even further by adding a while loop to it, that reruns the boxes game until the user enters “-1”.
For more details on programming in R, check out Nathaniel’s YaRrr! The Pirate’s Guide to R YaRrr! link
Hadley Wickham, covers programming in R and purrr in his book R for data science: R for Data Science link