Chapter 2 Higher Order Functions

One of the things we should demand from a powerful programming language is the ability to build abstractions by assigning names to common patterns and then to work in terms of the names directly.

Functions provide this ability.

Also, there are common programming patterns that recur in code, but are used with a number of different functions. These patterns can also be abstracted, by giving them names.

To express certain general patterns as named concepts, we will need to construct functions that can accept other functions as arguments or return functions as values.

Definition 2.1 Functions that manipulate functions are called higher-order functions.

2.1 Functions as arguments

Suppose we have the following R functions:

Sum of natural numbers up to n \[ \sum_{k=1}^n k \]

sum_naturals <- function(n){
    total <- 0
    for(k in 1:n){
        total <- total + k
    }
    return(total)
}

Sum of the cubes of natural numbers up to n.

\[ \sum_{k=1}^n k^3 \]

sum_cubes <- function(n){
    total <- 0
    for(k in 1:n){
        total <- total + k*k*k
    }
    return(total)
}

Sum of the series \(\frac{8}{1\times 3} + \frac{8}{5\times 7} + \frac{8}{9 \times 11}\cdots\) which converges to \(\pi\)

\[ \sum_{k=1}^n \frac{8}{(4k-3)\times(4k-1)} \]

pi_sum <- function(n){
    total <- 0
    for(k in 1:n){
        total <- total + 8/((4*k-3) * (4*k-1))
    }
    return(total)
}

These functions share a common underlying pattern. They are the same in most parts, only differing in:

function name
the function <func> applied on k inside the loop

<name> <- function(n){

  total <- 0
  for(k in 1:n){
    total <- total + <func>(k)
  }
  
  return(total)
  
}

If that is the case, we can use a function as an argument input of another function.

Suppose the following summation function, with an additional argument func.

summation <- function(n, func){
    total <- 0
    for(k in 1:n){
        total <- total + func(k)
    }
    return(total)
}

This implements the following summation:

\[ \sum_{k=1}^n f(k) \]

Examples

Summation of cube

\[ f(k) = k \times k \times k \]

cube <- function(k){
    k*k*k
}
summation(10, cube)

## [1] 3025

Approximation of \(\pi\)

\[ f(k) = \frac{8}{(4k-3)\times(4k-1)} \]

pi_term <- function(k){
    8/((4*k-3)*(4*k-1))
}

summation(1000000, pi_term)

## [1] 3.141592

2.2 Functions as return values

Another power of functions is that they can also return another function.

Function composition is the process of combining two or more functions to create a new function. The composed function applies one function to the result of another, enabling modular and reusable code.

That is, given two functions f(x) and g(x), we might want to define h(x)=f(g(x))

We can define function composition using our existing tools

Notice in the following function compose1, another function is defined instead, which is h.

compose1 <- function(f, g){
  h <- function(x){
    return(f(g(x)))
  }
  return(h)
}

For example, the goal is to compute the square of the input value + 1.

That is: \[ (x + 1)^2 \]

Let’s apply this function.

square <- function(x){
  return(x*x)
}

successor <- function(x){
  return(x + 1)
}

compose1 <- function(f, g){
  h <- function(x){
    return(f(g(x)))
  }
  return(h)
}

f <- function(x){
  # Never called.
  return(-x)
}

square_successor <- compose1(square, successor)

The square_successor will not have a return value yet, since it became a function as well.

square_successor

## function (x) 
## {
##     return(f(g(x)))
## }
## <environment: 0x0000028f7e739a80>

Now, applying this new defined function square_successor:

result <- square_successor(12)
result

## [1] 169

2.3 Functions as general methods

With higher-order functions, we begin to see a more powerful kind of abstraction: some functions express general methods of computation, independent of the particular functions they call.

Despite this conceptual extension of what a function means, our environment model of how to evaluate a call expression extends gracefully to the case of higher-order functions, without change.

When a user-defined function is applied to some arguments, the formal parameters are bound to the values of those arguments (which may be functions) in a new local frame.

`improve`: An iterative improvement algorithm

An iterative improvement algorithm begins with a “guess” of a solution to an equation.

It repeatedly applies an update function to improve that guess.

Apply a close comparison to check whether the current guess is “close enough” to be considered correct.

Once the guess is close enough, return it.

The following is an R function that implements this algorithm:

improve <- function(update, close, guess = 1){
    while(!close(guess)){
       guess <- update(guess)
    }
    return(guess)
}

The improve function is a general refinement loop, leaving the specific problem to the update and close functions.

Example: Golden Ratio

Consider the following example, which implements a general method for iterative improvement and uses it to compute (or estimate) the golden ratio.

The golden ratio \(\phi\) , often called “phi”, is a number near 1.6 that appears frequently in nature, art, and architecture.

Among the well-known properties of the golden ratio are that it can be computed by repeatedly summing the inverse of any positive number with 1, and that it is one less than its square.

\[\begin{align} \text{algebraic form}\\ \varphi &=\frac{1+\sqrt{5}}{2}\\ \text{property 1}\\ \varphi &= 1 + \frac{1}{1 + \frac{1}{1 + \ddots}}\\ &= 1+\frac{1}{\varphi} \\ \text{property 2}\\ \varphi &= \varphi^2-1 \end{align}\]

These two properties may be represented by the following functions to be used with improve:

golden_update <- function(guess){
  return(1 + 1/guess)
}

square_close_to_successor <- function(guess){
  return(approx_eq(guess, guess^2 - 1))
}

Above, we introduced the function approx_eq that is used to return TRUE if the two values are approximately equal to each other at some tolerance value.

approx_eq <- function(x, y, tolerance=1e-15){
  return(abs(x - y) < tolerance)
}

Calling improve with the arguments golden_update and square_close_to_successor will compute a finite approximation to the golden ratio.

improve(golden_update, square_close_to_successor)

## [1] 1.618034

Tracing the steps…

A local frame for improve is constructed with bindings for update, close, and guess.
In the body of improve, the name close is bound to square_close_to_successor, which is called on the initial value of guess.

As always, our new general method improve needs a test to check its correctness.

The golden ratio can provide such a test, because it also has an exact closed-form solution \(\varphi=(1+\sqrt{5})/2\), which we can compare to the iterative result of improve

phi <- 1/2 + sqrt(5)/2

approx_phi <- improve(golden_update, square_close_to_successor)

approx_eq(phi, approx_phi)

## [1] TRUE

2.4 Nested Definitions

Problems we might encounter with what we did earlier:

Global frame becomes cluttered with names of small functions, which must all be unique.
We are constrained by the fact that the update argument to improve must take exactly one argument.

Nested function definitions address both of these problems but require us to enrich our environment model.

Example: Square root of a number

Let’s consider a new problem: computing the square root of a number.

Repeated application of the following update converges to the square root of a number \(S\).

(Babylonian Method/Heron’s Method):

Begin with an arbitrary starting value \(x_0\) (close to the actual value of square root of \(S\), the better).
Let \(x_{n+1}\) be the average of \(x_n\) and \(S/x_n\)
Repeat step 2 until the desired accuracy is achieved.

It can also be represented as:

\[ \begin{align} x_0 &\approx \sqrt{S} \\ x_{n+1} &=\frac{1}{2}\left( x_n +\frac{S}{x_n}\right)\\ \sqrt{S} &=\lim_{n\rightarrow\infty}x_n \end{align} \]

Using the repeat loop, we can use this algorithm to estimate the square root of \(S=5\):

S <- 5

# initial guess for sqrt(5)
x <- 2
repeat{
    # update value of square root
    x <- (x + S/x)/2
    # exit if x^2 and S are close
    if(abs(x*x-S) <1e-15){
        break
    }
}

But the goal of this section is not to implement loops, but to use higher order functions, specifically nested definitions.

The repeated application of the following functions converges to the square root of \(S\).

average <- function(x, y){
  return((x + y)/2)
}

sqrt_update <- function(x, S){
  return(average(x, S/x))
}

This two-argument update function is incompatible with improve (they both take two arguments, not one), and it provides only a single update, while we really care about taking square roots by repeated updates.

The solution to both of these issues is to place function definitions inside the body of other definitions.

sqrt_cust <- function(S, guess){
    # first function defined inside a function definition
    sqrt_update <- function(x){
    return(average(x, S/x))
    }
    
    # second function defined inside a function definition
    sqrt_close <- function(x){
    return(approx_eq(x*x, S))
    }
    
    
    return(improve(sqrt_update, sqrt_close, guess))
}

Note here that the sqrt_update and sqrt_close will not be defined in the global environment.

sqrt_cust(5, guess = 2)

## [1] 2.236068

Lexical Scope

Locally defined functions also have access to the name bindings in the scope where they are defined.

This discipline of sharing names among nested definitions is called lexical scoping. Critically, the inner functions have access to the names in the environment where they are defined (not where they are called).

The scoping rules of a language determine how a value is associated with a free variable in a function.

Consider the function below:

func <- function(a, b){
    x + (a*b)
}

Free variables are not formal arguments and are not local variables (assigned inside the function body). In this example, x is a free variable.

R uses lexical scoping or static scoping.

An alternative to lexical scoping is dynamic scoping which is implemented by some languages.

We require two extensions to our environment model to enable lexical scoping.

Each user-defined function has a parent environment: the environment in which it was defined.
When a user-defined function is called, its local frame extends its parent environment.

To further illustrate, we have another example

y <- 10 
f <- function(x){
    y <- 2
    y^2 + g(x)
}

g <- function(x){
    x*y
}

The y in f()

The y in the definition of f() is not a formal parameter, but defined inside the body. This is a local variable.

The y in g()

On the other hand, in the new function g(), y is still not a formal argument and also not locally defined. This is a free variable.

Its value will be based on the parent environment where it was defined, not on the local frame of g().

Questions:

What is the value of f(3)?
What is the value of g(3)?

Practice Exercises

Applying a summary function

In R, simple summary statistics can be obtained using built-in functions such as sum, mean, sd, and median. Write a simple function, apply_function, that takes in two arguments: f, the function to be applied, and data which is a vector that contains the values to which the function will be applied to. To test your function, assign the following values to x.
```
x <- c(1, 2, -3, 6, 9)
```
Your function apply_function should return the value returned by the function f when applied to the vector data.
Applying a summary function on a trimmed vector

Modify your function apply_function such that it trims (removes) the maximum and minimum values in the vector data. Use the functions min and max and Boolean logic to accomplish this. Name your function apply_trimmed.

To exclude a value in a vector, you may use indexing, such as the following line of code:
```
# remove 2 from the vector x
x[x != 2]
```
```
## [1]  1 -3  6  9
```
Applying a function twice

Write an R function called apply_twice that takes function (f) and a value x as arguments. This function should apply f to x twice in succession and return the result.

For example, the expected output of apply_twice(square, 2) is 16, which performs \((2^2)^2 = 16\)
Generating a sequence

Write a higher-order function in R called generate_sequence that takes a function (f), a starting value (x), and a number of steps (n).

The function should generate a sequence by repeatedly applying f to the previous result, starting from x, and store all intermediate results in a vector.
1. Define the higher-order function generate_sequence(f, x, n).
2. The function should create a vector of length n+1 where
  - The first element is x
  - Each subsequent element is obtained by applying f to the previous element.
  - To create a vector, use the function c().
3. Return the generated sequence.
4. Test the function with the different sequences:
- Adding 3 (arithmetic progression)
```
add_three <- function(x) {x + 3}
generate_sequence(add_three, 0, 5)
# Expected Output: c(0, 3, 6, 9, 12, 15)
```
- Doubling numbers (geometric progression)
```
double <- function(x) {x * 2}
generate_sequence(double, 1, 5)
# Expected Output: c(1, 2, 4, 8, 16, 32)
```
Solving for x

Suppose you have the following equation:

\[ x = cos(x) \]

Apply the improve algorithm to solve for \(x\). Properly define the update and close functions prior to using them as arguments in the improve function.