5  Your first function

sqrt(9), c(1, 2, 3), sum(x > 70): every piece of R code so far has relied on functions that someone else wrote. This chapter shows you how to write your own.

5.1 Anatomy of a function

Here is a function that doubles its input:

double <- function(x) x * 2

Three parts. function(x) declares one parameter called x; the body, x * 2, is the expression R evaluates when you call the function; and double is the name you are binding this function to, the same way x <- 3 binds a name to a number. What you pass in when you call it is the argument.

Call it:

double(5)
#> [1] 10
double(c(1, 2, 3))
#> [1] 2 4 6

The second call works because * is vectorized (Section 4.4): the body x * 2 already operates element by element, so double handles vectors without a loop, without any extra effort on your part. But what happens when the body grows beyond a single expression?

A function with a longer body uses curly braces:

describe <- function(x) {
  n <- length(x)
  m <- mean(x)
  paste("n =", n, ", mean =", round(m, 2))
}

describe(c(72, 85, 61, 90, 78))
#> [1] "n = 5 , mean = 77.2"

The body is a block of expressions that R evaluates in order, returning the value of the last one. This is the same rule from Section 3.3: a block in curly braces is itself an expression, and its value is whatever the final line produces.

Exercises

  1. Write a function square that takes a number and returns its square.
  2. Write a function to_celsius that converts a Fahrenheit temperature to Celsius (formula: (f - 32) * 5/9). Test it with to_celsius(212).
  3. Does to_celsius(c(32, 72, 212)) work? Why?

5.2 Arguments and defaults

Functions can have multiple arguments, and some can have default values:

greet <- function(name, greeting = "Hello") {
  paste(greeting, name)
}

name is required; greeting is optional, defaulting to "Hello" if you leave it out:

greet("R")
#> [1] "Hello R"
greet("R", "Hey")
#> [1] "Hey R"

You can also pass arguments by name, in any order:

greet(greeting = "Bonjour", name = "R")
#> [1] "Bonjour R"

This is the same positional-versus-named distinction from Section 3.4, when we used log(8, base = 2). The rule is uniform across R because function calls all work the same way, whether you wrote the function or someone at R Core did.

Some functions accept a variable number of arguments using ... (pronounced “dot-dot-dot”):

shout <- function(...) {
  words <- c(...)
  toupper(paste(words, collapse = " "))
}

shout("hello", "world")
#> [1] "HELLO WORLD"

... collects every argument you pass, however many there are, and c(...) combines them into a vector. c(), paste(), and sum() all use ..., which is why c(1, 2, 3) and c(1, 2, 3, 4, 5, 6, 7) both work. More on ... in Chapter 23.

In the lambda calculus, a function takes exactly one argument. If you want two, you write a function that returns another function:

λx. λy. x + y

This reads: “a function that takes x and returns another function that takes y and returns x + y.” To add 3 and 5, you apply it twice: (λx. λy. x + y)(3)(5) = (λy. 3 + y)(5) = 3 + 5 = 8. Each application replaces one argument, one step at a time. This is called currying, and R can do it too (Chapter 20). ... is the practical alternative: a single function accepts any number of arguments and decides what to do with them inside.

Exercises

  1. Write a function power with arguments x and n = 2 that computes x^n. Test it with power(3) and power(3, 3).
  2. What happens if you call greet() with no arguments? Read the error message.
  3. Rewrite power so the default exponent is 3. What does power(2) return now?

5.3 Return values

R returns the value of the last expression in the body. No keyword needed:

add_one <- function(x) x + 1
add_one(10)
#> [1] 11

The body is x + 1. R evaluates it and returns the result. A function call is an expression that evaluates to a value, like every other expression in R (Section 3.3).

When you call add_one(10), R substitutes 10 for x in the body, turning x + 1 into 10 + 1, then evaluates it. The lambda calculus calls this beta reduction: (λx. x + 1)(10) reduces to 10 + 1. Every function call in R is a beta reduction.

R does have a return() function, and it is useful for early exits:

safe_log <- function(x) {
  if (any(x <= 0)) return(NA)
  log(x)
}

safe_log(5)
#> [1] 1.609438
safe_log(-1)
#> [1] NA

If x has any non-positive values, the function exits immediately and returns NA. Otherwise, it falls through to log(x).

TipOpinion

Use return() only for early exits like safe_log above. If the function runs to the end, let the last expression be the return value; writing return(x) on the last line adds nothing but noise.

One subtlety: assignment returns its value invisibly. When a function’s last line is an assignment, R returns the value but does not print it:

quiet <- function(x) {
  result <- x + 1
}

quiet(5)

Nothing printed. The function returned 6, but invisibly. If you want to see it:

(quiet(5))
#> [1] 6

Or store it:

y <- quiet(5)
y
#> [1] 6

The general rule: if you want your function to show a result, make sure the last expression is the value itself, not an assignment. But there is a deeper question lurking here: if a function is just a value, can you treat it the way you treat a number or a string?

Exercises

  1. What does this function return? Predict, then check:

    mystery <- function(x) {
      a <- x + 1
      b <- a * 2
      b
    }
    mystery(3)
  2. What about this one?

    mystery2 <- function(x) {
      a <- x + 1
      b <- a * 2
    }
    mystery2(3)
  3. Write a function clamp that takes a value x, a lo, and a hi, and returns lo if x < lo, hi if x > hi, and x otherwise. Use return() for the early exits.

5.4 Functions are just values

Type sqrt without parentheses and R prints the function itself, not a result. A function is a value, the same kind of thing as a number or a string: you can assign it, put it in a list, or interrogate its type.

f <- double
f(5)
#> [1] 10
is.function(double)
#> [1] TRUE
typeof(double)
#> [1] "closure"

double is a name bound to a value, and that value happens to be a function. f <- double gives the same function a second name, which is exactly what happened in Section 3.2 when y <- x gave the number 3 a second name.

You can even create a function without giving it a name at all:

(function(x) x + 1)(10)
#> [1] 11

function(x) x + 1 creates a function; the parentheses around it and (10) immediately call it with the argument 10. This is an anonymous function, λx. x + 1 applied to 10, and Chapter 19 uses them heavily.

You can do anything with a function that you can do with a number: store it in a variable, pass it to another function, return it as a result. This is what McCarthy built into Lisp in 1958 (Section 2.1), and R inherits it directly. So functions are values, and values live somewhere. Where do the variables inside a function live?

Exercises

  1. Assign sqrt to a new name s. Does s(16) work?
  2. What does is.function(42) return? What about is.function(is.function)?
  3. Run (function(a, b) a + b)(3, 4). What happens?

5.5 Scope: a first look

What happens to variables created inside a function?

add_ten <- function(x) {
  y <- 10
  x + y
}

add_ten(5)
#> [1] 15

y was created inside add_ten. Can you access it from outside?

y
#> [1] 6

No. Variables created inside a function live in the function’s own local environment and vanish when the function finishes. This is what makes functions safe: they do not leak their internal state into your workspace. But the boundary is not symmetric.

multiplier <- 3

scale <- function(x) {
  x * multiplier
}

scale(10)
#> [1] 30

scale found multiplier in the workspace because it was not defined inside the function body. R looks for names first in the local environment, then in the environment where the function was defined, a strategy called lexical scoping that R inherited from Scheme (Section 2.2).

If a function can read outer variables, can it change them?

n <- 0

increment <- function() {
  n <- n + 1
  n
}

increment()
#> [1] 1
n
#> [1] 0

increment created its own local n, set it to 1, and returned it, but the n in the workspace is still 0. The function read the outer n without modifying it; it created a new local binding with the same name. This is R’s copy-on-modify philosophy: values do not change, you create new values.

The practical rule is simple: what happens inside a function stays inside the function. The full machinery behind this (environments, closures, and the scope chain) is in Chapter 18.

Exercises

  1. Predict the output, then run it:

    x <- 100
    f <- function() {
      x <- 1
      x
    }
    f()
    x
  2. What does this return?

    x <- 10
    g <- function(x) {
      x + 1
    }
    g(20)

    Is it 11 or 21? Why?

  3. Can a function call another function you’ve defined? Try it: write triple <- function(x) x * 3, then write triple_and_add <- function(x) triple(x) + 1. Does triple_and_add(5) work?

5.6 Guarding your inputs

What happens when someone passes the wrong type to your function?

to_celsius <- function(f) (f - 32) * 5/9
to_celsius("hot")
#> Error in `f - 32`:
#> ! non-numeric argument to binary operator

The error says non-numeric argument to binary operator, which is technically accurate but unhelpful: you have to read the function body to figure out what went wrong. Better to catch the problem at the door:

to_celsius <- function(f) {
  stopifnot(is.numeric(f))
  (f - 32) * 5/9
}

to_celsius("hot")
#> Error in `to_celsius()`:
#> ! is.numeric(f) is not TRUE

stopifnot() checks a condition and stops with an error if it is FALSE. The message now tells you exactly what failed: is.numeric(f) is not TRUE. You can stack multiple checks in the same call:

safe_divide <- function(x, y) {
  stopifnot(is.numeric(x), is.numeric(y), y != 0)
  x / y
}

safe_divide(10, 0)
#> Error in `safe_divide()`:
#> ! y != 0 is not TRUE

A function that fails immediately with a clear message saves you from chasing wrong outputs later. Chapter 25 covers input validation in depth; for now, stopifnot() at the top of your functions is enough. But validation only protects against bad inputs. What about the structure of the computation itself?

Exercises

  1. Add a stopifnot() to your power function that checks x is numeric and n is numeric.
  2. What happens if you write stopifnot(is.numeric(x), length(x) == 1) and pass a vector? Try it.

5.7 Recursion

A function can call other functions. Can it call itself?

factorial <- function(n) {
  if (n == 0) return(1)
  n * factorial(n - 1)
}

factorial(5)
#> [1] 120

factorial(5) calls factorial(4), which calls factorial(3), and so on, until factorial(0) returns 1 and the whole chain collapses back:

factorial(5)
= 5 * factorial(4)
= 5 * 4 * factorial(3)
= 5 * 4 * 3 * factorial(2)
= 5 * 4 * 3 * 2 * factorial(1)
= 5 * 4 * 3 * 2 * 1 * factorial(0)
= 5 * 4 * 3 * 2 * 1 * 1
= 120

This is recursion: a function defined in terms of itself. Look at the unfolding carefully. Each line replaces a function call with its value, the same term replacement from Section 4.2, and no variable is being updated, no counter is ticking. The answer was implicit in the original expression factorial(5); evaluation simply made it explicit.

Figure 5.1: Recursive call stack for factorial(5): each frame waits for the one below it. No variable changes; each frame holds its own n.

This is how the lambda calculus handles repetition. There is no loop counter; the function calls itself with a smaller argument until it stops.

In lambda calculus notation, factorial looks like this:

fact = λn. if n = 0 then 1 else n * fact(n - 1)

The self-reference is informal shorthand; pure lambda calculus has no names, so fact cannot refer to itself directly, and Section 30.5 shows how recursion works without them. The if/then/else is shorthand too. The lambda calculus has only variables, functions, and application, but Church showed that true, false, and branching can be encoded as functions that choose between their arguments:

TRUE  = λx. λy. x      # a function that takes two arguments and returns the first
FALSE = λx. λy. y      # a function that takes two arguments and returns the second
IF    = λcond. λthen. λelse. cond(then)(else)

IF hands both options to the condition and lets it choose, so IF(TRUE)(1)(2) becomes TRUE(1)(2) becomes (λx. λy. x)(1)(2) becomes 1. Section 30.2 builds these in R.

With the shorthand in place, the reduction works as before, one substitution at a time:

fact(3)                                     # apply fact with n = 3
= (λn. if n=0 then 1 else n * fact(n-1))(3) # replace fact with its definition
= if 3=0 then 1 else 3 * fact(3-1)          # substitute n = 3
= 3 * fact(2)                               # 3≠0, so take the else branch
= 3 * 2 * fact(1)                           # same pattern: n=2, substitute, reduce
= 3 * 2 * 1 * fact(0)                       # n=1
= 3 * 2 * 1 * 1                             # n=0, take the then branch: 1
= 6

R’s function(n) is the same as λn, and this reduction is what R does when you call factorial(3): substitute 3 for n, evaluate the body, hit factorial(2), substitute again, keep going until the base case.

fibonacci <- function(n) {
  if (n <= 1) return(n)
  fibonacci(n - 1) + fibonacci(n - 2)
}

fibonacci(10)
#> [1] 55

fibonacci(30) is already slow, because the function recomputes the same values over and over: fibonacci(10) calls fibonacci(9) and fibonacci(8), but fibonacci(9) also calls fibonacci(8), so that subtree gets computed twice. Further down, fibonacci(7) gets computed three times, fibonacci(6) five times, and the work grows exponentially. Is there a way to keep the recursive structure without paying the exponential cost?

A smarter version avoids this by carrying the last two values forward as arguments:

fibonacci_fast <- function(n, a = 0, b = 1) {
  if (n == 0) return(a)
  fibonacci_fast(n - 1, b, a + b)
}

fibonacci_fast(10)
#> [1] 55
fibonacci_fast(50)
#> [1] 12586269025

The naive version cannot compute fibonacci(50) in any reasonable time; this one handles it instantly. Instead of branching into two recursive calls, fibonacci_fast makes one call per step, passing the accumulated values along. Watch the reduction:

fibonacci_fast(4, 0, 1)
= fibonacci_fast(3, 1, 1)    # n=4→3, a=0→1, b=1→0+1=1
= fibonacci_fast(2, 1, 2)    # n=3→2, a=1→1, b=1→1+1=2
= fibonacci_fast(1, 2, 3)    # n=2→1, a=1→2, b=2→1+2=3
= fibonacci_fast(0, 3, 5)    # n=1→0, a=2→3, b=3→2+3=5
= 3                           # n=0, return a

No variable is being mutated; each call just substitutes new values for n, a, and b and hands them to the next call. This form is called tail recursion because the recursive call is the last thing the function does.

Any function with multiple recursive calls can be rewritten as a single recursive call with extra arguments that carry the intermediate results. For Fibonacci, the state to pass forward is the last two numbers; for a tree traversal, it might be an accumulator list.

R has a limit on how deep recursion can go (around 1000 calls by default), and for most day-to-day R work, vectorized operations and functions like map() (Chapter 19) are better tools.

But there is something strange about factorial that has nothing to do with performance: its base case and inductive step look exactly like a proof by mathematical induction.

Exercises

  1. Trace through factorial(3) by hand, writing out each step like the example above.
  2. Write a recursive function sum_to that adds up all integers from 1 to n. Check: sum_to(100) should return 5050.
  3. What happens if you call factorial(-1)? Why? Add a stopifnot() to guard against it.

5.8 Programs and proofs

Look at factorial again:

factorial <- function(n) {
  if (n == 0) return(1)
  n * factorial(n - 1)
}

Two parts: a base case (when n = 0, return 1) and an inductive step (assume factorial(n - 1) gives the right answer, and use it to compute factorial(n)). If you have seen mathematical induction, you recognize the shape immediately. A proof by induction has a base case and an inductive step too. The resemblance is not a coincidence.

In the 1930s, a logician working on how functions combine noticed something strange: the rules for combining types looked identical to the rules for combining propositions in formal logic. A function that takes an A and returns a B had the same structure as a proof that “if A, then B.” His name was Haskell Curry, and the Haskell programming language is named after him. Three decades later William Howard made the observation precise: every type corresponds to a logical proposition, and every program of that type corresponds to a proof of that proposition. A recursive function that terminates is a proof by induction. The result is called the Curry-Howard correspondence.

Writing factorial constructs a proof by induction that the product 1 * 2 * ... * n exists for every non-negative integer. The transformation from double-recursive fibonacci to tail-recursive fibonacci_fast is the same move a logician makes when they strengthen an inductive hypothesis to make a proof go through: the extra arguments (a and b) are the stronger hypothesis. Double induction reduces to simple induction with stronger invariants, which is why the same trick turns a slow recursive function into a fast one.

Church, Curry, and Turing were all asking, before any computer existed, what could be proved, and the answer turned out to be identical to what could be computed.

Look at where the proof lives in factorial. The definition says that factorial(n) and n * factorial(n - 1) are the same expression for arbitrary n, an equivalence that holds by term replacement alone, without evaluating either side. It is true by definition. The base case says factorial(0) is 1, also by definition. Those two statements are the entire proof: you show that the program in state n and the program in state n - 1 are the same program up to one substitution, and that the chain has a starting point. That is a proof by induction, and it covers every natural number at once. In lambda calculus notation:

fact = λn. if n = 0 then 1 else n * fact(n - 1)

fact is the function. fact(n) is the function applied to some n. Replacing the name fact with its definition gives:

fact(n)
= (λn. if n = 0 then 1 else n * fact(n - 1))(n)   # replace fact with its definition
= if n = 0 then 1 else n * fact(n - 1)              # substitute n

You might think: of course they are the same, we just wrote the definition. But fact and fact(n) are not the same thing. fact is a definition, a rule that says what to do with an argument; fact(n) is that rule applied to a specific (but arbitrary) n.

The reduction above shows what happens when you apply the rule: you get a term that contains fact(n - 1), the same program applied to a smaller argument. That relationship between fact(n) and fact(n - 1) is not something we assumed or proved; it fell out of one substitution. A tautology. And yet that tautology, together with the base case, constitutes a proof by induction that fact produces the right answer for every natural number. When n = 0, the expression reduces to 1. The substitution gives you the step from n - 1 to n.

This is what the Curry-Howard correspondence means in practice: the equivalence of two program states is a tautology that requires only substitution, and a chain of such tautologies with a starting point is a proof.

The correspondence is the foundation of formal verification: software that proves other software is correct. The first practical tool built on it was Coq, developed at INRIA in the late 1980s. In Coq, every program must come with a machine-checkable proof that it does what its type signature claims; if the proof has a gap, the code does not compile. This makes Coq painful for ordinary programming but powerful for building software where correctness matters more than convenience.

In 2013, Leonardo de Moura at Microsoft Research started Lean, a newer proof assistant with a more accessible syntax and a growing library of formalized mathematics called Mathlib. In Lean, a proof that adding zero to any natural number gives back the same number looks like this:

theorem add_zero (n : Nat) : n + 0 = n := by
  induction n with
  | zero => rfl
  | succ n ih => simp [Nat.add_succ, ih]

The structure is identical to factorial: a base case (zero), an inductive step (succ n), and the result from the smaller case (ih, the induction hypothesis). If the proof is wrong, the code does not compile. The proof is a program, and the program is a proof.

The seL4 microkernel (used in military helicopters and autonomous vehicles) has a machine-checked proof that certain classes of errors in its code cannot occur, for any input. Not “we haven’t found a bug yet”: a mathematical guarantee, verified by a computer. The CompCert C compiler has a similar proof that it will never miscompile your code. Both proofs were written as programs in Coq. Traditional testing covers a finite set of inputs; a proof by induction, because it operates on the structure of the program rather than on specific values, covers all of them.

R checks none of this; nothing stops you from writing a recursive function with no base case. But the structure is there in every recursive function you write, base case and inductive step, the same shape as a proof that holds for all natural numbers.