Here’s a function with a single argument. The function simply adds one to the argument. Putting out on the last line will make this the function’s output.
addone <- function(x){
out <- x+1
out
}Let’s try this function out with a scalar and then a vector.
addone(3)## [1] 4these <- 1:10
addone(these)## [1] 2 3 4 5 6 7 8 9 10 11If we want, we can simplify by just calculating the value (and not saving it with the arrow).
addone <- function(x){
x+1
}We can now test it out, for fun.
addone(3)## [1] 4these <- 1:10
addone(these)## [1] 2 3 4 5 6 7 8 9 10 11Our output doesn’t have to be a scalar. Here’s an example of a function that creates and outputs a diagonal matrix with the number of rows specified by m.
makemat <- function(m){
diag(m)
}
makemat(3)## [,1] [,2] [,3]
## [1,] 1 0 0
## [2,] 0 1 0
## [3,] 0 0 1Now, YOU write a function that evaluates 2x+1 for any given value of x.
Oftentimes, you will want to output more than a single object. A nice way to output several objects is by employing a list. Do all your calculations, then create a list using the list() function. I also recommend naming the objects in your last, as I show below.
domath <- function(x){
expit <- exp(x)
doubleit <- 2*x
out <- list(expit = expit, doubleit = doubleit)
out
}Test it out!
domath(2)## $expit
## [1] 7.389056
##
## $doubleit
## [1] 4Save the output and do stuff with it!
out <- domath(2)
out## $expit
## [1] 7.389056
##
## $doubleit
## [1] 4If you’ve said it, you can dig out the pieces of output in a few ways. First, you could use indexing (square brackets). This gives the first item in the list:
out[1]## $expit
## [1] 7.389056and here’s the second:
out[2]## $doubleit
## [1] 4A second way you can dig out the items of the list is by using the names. First, figure out what the names are.
names(out)## [1] "expit" "doubleit"Then use them with the dollar sign operator.
out$expit## [1] 7.389056out$doubleit## [1] 4Now, YOU add a third computation to the list in the domath() function.
Want to use more than one argument? Of course you do. Don’t be stupid. Let’s solve a quadratic equation for fun. First, let’s define our quadratic function so we can be clear on what a, b and c are: \[ax^2 + bx +c =0\].
solveit <- function(a, b, sea){
top <- -b +c(-1,1) * sqrt(b^2 - 4*a*sea)
bottom <- 2*a
x <- top/bottom
x
}Now let’s test it first with a baby example. When does \(4x^2 =0\)? For \(x = 0\). Let’s see if it gives the right answer.
solveit(a=4, b=0, sea=0)## [1] 0 0When does \(10x^2 - 9x -6 =0\)?
solveit(a=10, b = -9, sea = -6)## [1] -0.4458236 1.3458236If we check this, it’s right. (Phew.)
top <- 9 +c(-1,1) * sqrt(321)
bottom <- 2*10
top/bottom## [1] -0.4458236 1.3458236Did you notice that we named the arguments in the function when we called the function? It’s nice to do for multi-argument functions, especially so we can enter the arguments out of order. For example, compare the following output to the previous output.
solveit(a=10, sea = -6, b = -9)## [1] -0.4458236 1.3458236If we take away the labels, we get totally different answers.
solveit(10, -6, -9)## [1] -0.6949874 1.2949874I was reading about carpal tunnel prevention and learned that BMI is a risk factor. Naturally, I was curious about my BMI and calculated it using the NHLBI’s BMI calculator. I don’t like black boxes though, so I wanted to find the equation, which of course led to me deciding I should program it so you can see a non-trivial function in action. BMI is calculated using height and weight, but in meters and kilograms. So that you don’t have to translate your height and weight to meters and kilograms, the following function takes the weight and height in pounds, feet, and inches. I could have just outputted the resulting BMI, but I thought I’d also output the height and weight in meters and kilograms to check my work (and educate myself on my height and weight in these units!).
bmi <- function(pounds, feet, inches) {
weight <- pounds * .453 #calculates weight in kg
height <- (30.48 * feet + 2.54 * inches)/100 #calculates height in m
out <- weight/(height^2) #this is your bmi
list(weight = weight, height = height, result = out)
}Testing it out for a person who is 5’9" and 136 pounds produces a BMI of 20.06, as shown below.
bmi(136, 5, 9)## $weight
## [1] 61.608
##
## $height
## [1] 1.7526
##
## $result
## [1] 20.05725We had a function do many computations and then output them all. What if we want a function that holds many functions so that we can potentially execute just a subset of the functions? I do this in funlist below.
funlist <-function()
{
square <- function(x) x^2 #squares your term
addone <- function(x) x+1 #adds one to your term
lessbaby <- function(x){ #adds one then squares the result
temp <- funlist()$addone(x)
funlist()$square(temp)
}
out<-list(square = square, addone = addone, lessbaby = lessbaby)
out
}Here’s a way to use each function (one step process). The following runs each function in the function list using 2 as its argument.
funlist()$square(2)## [1] 4funlist()$addone(2)## [1] 3funlist()$lessbaby(2)## [1] 9If you want to refer to functions using indexing rather than function name, you can do that! The following grabs first function, applies it to the number 2. Woof; that’s ugly. Someone find a way to clean that up!
funlist()[[1]](2) ## [1] 4Does it seem silly? To have a function contain a list of functions? It might be silly, but it’s used in the very popular glm() and in the less-popular-but-still-awesome glmm().
If you need to execute a chunk of code iteratively, a loop might be a good idea.
Say we want to collect a roll a fair six-sided die m times and record the result. We could use a for loop.
m <- 100 # number of times you roll the die
roll <- rep(0,m)
for(i in 1:m){
roll[i] <- sample(1:6, size = 1)
}Let’s roll a die m times and calculate the proportion of each side.
table(roll)/m## roll
## 1 2 3 4 5 6
## 0.23 0.11 0.19 0.17 0.19 0.11Note: this is an illustrative example to demonstrate for loops in an easy-to-understand way. If you actually want m die rolls, you can do this more quickly by sampling m times with replacement:
results <- sample(1:6, size = m, replace = TRUE)
table(results)/m## results
## 1 2 3 4 5 6
## 0.18 0.22 0.18 0.14 0.14 0.14While loops will keep iterating the loop while something is true. Let’s keep rolling til we get a 6, and record the result the whole time.
rolls <- c() #bad practice, sorry!
nosix <- TRUE # here's an indicator
i <- 1
while(nosix == TRUE){
rolls[i] <- sample(1:6, size = 1)
if(rolls[i] == 6){ nosix <- FALSE}
i <- i + 1
}
rolls## [1] 3 4 3 3 1 1 2 4 1 1 5 5 2 2 4 5 6