Chapter 7
Algorithms and Numerical Methods

Although the percentage of time that a quant spends on programming varies with the job function (e.g., quant analyst/researcher versus quant developer) and firm culture, a typical quant generally devotes part of his or her time to implementing models through programming. Therefore, programming skill test is often an inherent part of the quantitative interview.

To a great extent, the programming problems asked in quantitative interviews are similar to those asked in technology interviews. Not surprisingly, many of these problems are platform- or language-specific. Although C++ and Java still dominate the market, we’ve seen a growing diversification to other programming languages such as Matlab, SAS, S-Plus, and R. Since there are many existing books and websites dedicated to technology interviews, this chapter will not give a comprehensive review of programming problems. Instead, it discusses some algorithm problems and numerical methods that are favorite topics of quantitative interviews.

7.1. Algorithms

In programming, the analysis of algorithm complexity often uses asymptotic analysis that ignores machine-dependent constants and studies the running time $T(n)$—the number of primitive operations such as addition, multiplication, and comparison—as the number of inputs $n\to\infty$.¹

Three of the most important notations in algorithm complexity are big-$O$ notation, $\Omega$ notation and $\Theta$ notation:

$O(g(n))=\{f(n):$ there exist positive constants $c$ and $n_0$ such that $0\le f(n)\le cg(n)$ for all $n\ge n_0\}$. It is the asymptotic upper bound of $f(n)$.

$\Omega(g(n))=\{f(n):$ there exist positive constants $c$ and $n_0$ such that $0\le cg(n)\le f(n)$ for all $n\ge n_0\}$. It is the asymptotic lower bound of $f(n)$.

$\Theta(g(n))=\{f(n):$ there exist positive constants $c_1$, $c_2$, and $n_0$ such that $c_1g(n)\le f(n)\le c_2g(n)$ for all $n\ge n_0\}$. It is the asymptotic tight bound of $f(n)$.

Besides notations, it is also important to explain two concepts in algorithm complexity:

Worst-case running time $W(n)$: an upper bound on the running time for any $n$ inputs.

Average-case running time $A(n)$: the expected running time if the $n$ inputs are randomly selected.

For many algorithms, $W(n)$ and $A(n)$ have the same $O(g(n))$. But as we will discuss in some problems, they may well be different and their relative importance often depends on the specific problem at hand.

A problem with $n$ inputs can often be split into $a$ subproblems with $n/b$ inputs in each subproblem. This paradigm is commonly called divide-and-conquer. If it takes $f(n)$ primitive operations to divide the problem into subproblems and to merge the solutions of the subproblems, the running time can be expressed as a recurrence equation $T(n)=aT(n/b)+f(n)$, where $a\ge1,b>1$, and $f(n)\ge0$.

The master theorem is a valuable tool in finding the tight bound for recurrence equation $T(n)=aT(n/b)+f(n)$: If $f(n)=O(n^{\log_ba-\varepsilon})$ for some constant $\varepsilon>0$, $T(n)=\Theta(n^{\log_ba})$, since $f(n)$ grows slower than $n^{\log_ba}$. If $f(n)=\Theta(n^{\log_ba}\log^kn)$ for some $k\ge0$, $T(n)=\Theta(n^{\log_ba}\log^{k+1}n)$, since $f(n)$ and $n^{\log_ba}$ grow at similar rates. If $f(n)=\Omega(n^{\log_ba+\varepsilon})$ for some constant $\varepsilon>0$, and $af(n/b)\le cf(n)$ for some constant $c\lt1$, $T(n)=\Theta(f(n))$, since $f(n)$ grows faster than $n^{\log_ba}$.

Let’s use binary search to show the application of the master theorem. To find an element in an array, if the numbers in the array are sorted ($a_1\le a_2\le\cdots\le a_n$), we can use binary search: The algorithm starts with $a_{\lfloor n/2\rfloor}$. If $a_{\lfloor n/2\rfloor}=x$, the search stops. If $a_{\lfloor n/2\rfloor}>x$, we only need to search $a_1,\cdots,a_{\lfloor n/2-1\rfloor}$. If $a_{\lfloor n/2\rfloor}\lt x$, we only need to search $a_{\lfloor n/2+1\rfloor},\cdots,a_n$. Each time we can reduce the number of elements to search by half after making one comparison. So we have $a=1$, $b=2$, and $f(n)=1$. Hence, $f(n)=\Theta(n^{\log_21}\log^0n)$ and the binary search has complexity $\Theta(\log n)$.

Number swap

How do you swap two integers, $i$ and $j$, without using additional storage space?

Solution: Comparison and swap are the basic operations for many algorithms. The most common technique for swap uses a temporary variable, which unfortunately is forbidden in this problem since the temporary variable requires additional storage space. A simple 173 mathematic approach is to store the sum of $i$ and $j$ first, then extract $i$’s value and assign it to $j$ and finally assign $j$’s value to $i$. The implementation is shown in the following code:²

void swap(int &i, int &j) {
    i = i + j;  //store the sum of i and j
    j = i - j;  //change j to i's value
    i = i - j;  //change i to j's value
}

An alternative solution uses bitwise XOR (^) function by taking advantage of the fact that $x\mathbin{\wedge}x=0$ and $0\mathbin{\wedge}x=x$:

void swap(int &i, int &j){
    i = i ^ j;
    j = j ^ i; //j = i ^ (j ^ i) = i
    i = i ^ j; //i = (i ^ j) ^ i = j
}

Unique elements

If you are given a sorted array, can you write some code to extract the unique elements from the array? For example, if the array is [1, 1, 3, 3, 3, 5, 5, 5, 9, 9, 9, 9], the unique elements should be [1, 3, 5, 9].

Solution: Let $a$ be an $n$-element sorted array with elements $a_0\le a_1\le\cdots\le a_{n-1}$. Whenever we encounter a new element $a_i$ in the sorted array, its value is different from its previous element ($a_i\ne a_{i-1}$). Using this property we can easily extract the unique elements. One implementation in C++ is shown as the following function:³

template <class T> vector<T> unique(T a[], int n) {
    vector<T> vec; // vector used to avoid resizing problem
    vec.reserve(n); //reserver to avoid reallocation
    vec.push_back(a[0]);
    for(int i=1; i<n; ++i) {

        if(a[i] != a[i-1])
            vec.push_back(a[i]);
    }
    return vec;
}

Horner’s algorithm

Write an algorithm to compute $y=A_0+A_1x+A_2x^2+A_3x^3+\cdots+A_nx^n$.

Solution: A naïve approach calculates each component of the polynomial and adds them up, which takes $O(n^2)$ number of multiplications. We can use Horner’s algorithm to reduce the number of multiplications to $O(n)$. The algorithm expresses the original polynomial as $y=((((A_nx+A_{n-1})x+A_{n-2})x+\cdots+A_2)x+A_1)x+A_0$ and sequentially calculate $B_n=A_n,\ B_{n-1}=B_nx+A_{n-1},\cdots,B_0=B_1x+A_0$. We have $y=B_0$ with at most $n$ multiplications.

Moving average

Given a large array $A$ of length $m$, can you develop an efficient algorithm to build another array containing the $n$-element moving average of the original array ($B_1,\cdots,B_{n-1}=NA,\ B_i=(A_{i-n+1}+A_{i-n+2}+\cdots+A_i)/n,\ \forall i=n,\cdots,m$)?

Solution: When we calculate the moving average of the next $n$ consecutive numbers, we can reuse the previously computed moving average. Just multiply that average by $n$, subtract the first number in that moving average and then add the new number, and you have the new sum. Dividing the new sum by $n$ yields the new moving average. Here is the pseudo-code for calculating the moving average:

S = A[1] + ··· + A[n]; B[n] = S/n;
for (i=n+1 to m) { S = S - A[i-n] + A[i]; B[i] = S/n; }

Sorting algorithm

Could you explain three sorting algorithms to sort $n$ distinct values $A_1,\cdots,A_n$ and analyze the complexity of each algorithm?

Solution: Sorting is a fundamental process that is directly or indirectly implemented in many programs. So a variety of sorting algorithms have been developed for different 175 purposes. Here let’s discuss three such algorithms: insertion sort, merge sort and quick sort.

Insertion sort: Insertion sort uses an incremental approach. Assume that we have sorted subarray $A[1,\cdots,i-1]$. We insert element $A_i$ into the appropriate place in $A[1,\cdots,i-1]$, which yields sorted subarray $A[1,\cdots,i]$. Starting with $i=1$ and increases $i$ step by step to $n$, we will have a fully sorted array. For each step, the expected number of comparisons is $i/2$ and the worst-case number of comparisons is $i$. So we have

Merge sort: Merge sort uses the divide-and-conquer paradigm. It divides the array into two subarrays each with $n/2$ items and sorts each subarray. Unless the subarray is small enough (with no more than a few elements), the subarray is again divided for sorting. Finally, the sorted subarrays are merged to form a single sorted array.

The algorithm can be expressed as the following pseudocode:

mergesort(A, beginindex, endindex)
if beginindex < endindex
  then centerindex ← (beginindex + endindex)/2
       merge1 <- mergesort(A, beginindex, centerindex)
       merge2 <- mergesort(A, centerindex + 1, endindex)
       merge(merge1, merge2)

The merge of two sorted arrays with $n/2$ elements each into one array takes $\Theta(n)$ primitive operations. The running time $T(n)$ follows the following recursive function:

Applying the master theorem to $T(n)$ with $a=2$, $b=2$, and $f(n)=\Theta(n)$, we have $f(n)=\Theta(n^{\log_ba}\log^0n)$. So $T(n)=\Theta(n\log n)$. For merge sort, $A(n)$ and $W(n)$ are the same as $T(n)$.

Quicksort: Quicksort is another recursive sorting method. It chooses one of the elements, $A_i$, from the sequence and compares all other values with it. Those elements smaller than $A_i$ are put in a subarray to the left of $A_i$; those elements larger than $A_i$ are put in a subarray to the right of $A_i$. The algorithm is then repeated on both subarrays (and any subarrays from them) until all values are sorted.

In the worst case, quicksort requires the same number of comparisons as the insertion sort. For example, if we always choose the first element in the array (subarray) and compare all other elements with it, the worst case happens when $A_1,\cdots,A_n$ are already sorted. In such cases, one of the subarray is empty and the other has $n-1$ element. Each step only reduces the subarray size by one. Hence, $W(n)=\Theta(\sum_{i=1}^ni)=\Theta(n^2)$.

To estimate the average-case running time, let’s assume that the initial ordering is random so that each comparison is likely to be any pair of elements chosen from $A_1,\cdots,A_n$. If we suspect that the original sequence of elements has a certain pattern, we can always randomly permute the sequence first with complexity $\Theta(n)$ as explained in the next problem. Let $\widetilde A_p$ and $\widetilde A_q$ be the $p$th and $q$th element ($1\le p\lt q\le n$) in the final sorted array. There are $q-p+1$ numbers between $\widetilde A_p$ and $\widetilde A_q$. The probability that $\widetilde A_p$ and $\widetilde A_q$ is compared is the probability that $\widetilde A_q$ is compared with $\widetilde A_p$ before $\widetilde A_{p+1},\cdots,$ or $\widetilde A_{q-1}$ is compared with either $\widetilde A_p$ or $\widetilde A_q$ (otherwise, $\widetilde A_p$ and $\widetilde A_q$ are separated into different subarrays and will not be compared), which happens with probability $\displaystyle P(p,q)=\frac{2}{q-p+1}$ (you can again use the symmetry argument to derive this probability).

The total expected number of comparison is $\displaystyle A(n)=\sum_{q=2}^n\sum_{p=1}^{q-1}P(p,q)=\sum_{q=2}^n\sum_{p=1}^{q-1}\left(\frac{2}{q-p+1}\right)=\Theta(n\lg n)$.

Although theoretically quicksort can be slower than merge sort in the worst cases, it is often as fast as, if not faster than, merge sort.

Random permutation

A. If you have a random number generator that can generate random numbers from either discrete or continuous uniform distributions, how do you shuffle a deck of 52 cards so that every permutation is equally likely?

Solution: A simple algorithm to permute $n$ elements is random permutation by sorting. It assigns a random number to each card and then sorts the cards in order of their assigned random numbers.⁴ By symmetry, every possible order (out of $n!$ possible ordered sequences) is equally likely. The complexity is determined by the sorting step, so the 177 running time is $\Theta(n\log n)$. For a small $n$, such as $n=52$ in a deck of cards, the complexity $\Theta(n\log n)$ is acceptable. For large $n$, we may want to use a faster algorithm known as the Knuth shuffle. For $n$ elements $A[1],\cdots,A[n]$, the Knuth shuffle uses the following loop to generate a random permutation:

for (i=1 to n) swap(A[i], A[Random(i, n)]),

where Random($i,n$) is a random number from the discrete uniform distribution between $i$ and $n$.

The Knuth shuffle has a complexity of $\Theta(n)$ and an intuitive interpretation. In the first step, each of the $n$ cards has equal probability of being chosen as the first card since the card number is chosen from the discrete uniform distribution between 1 and $n$; in the second step, each of the remaining $n-1$ cards elements has equal probability of being chosen as the second card; and so on. So naturally each ordered sequence has $1/n!$ probability.

B. You have a file consisting of characters. The characters in the file can be read sequentially, but the length of the file is unknown. How do you pick a character so that every character in the file has equal probability of being chosen?

Solution: Let’s start with picking the first character. If there is a second character, we keep the first character with probability 1/2 and replace the pick with the second character with probability 1/2. If there is a third character, we keep the pick (from the first two characters) with probability 2/3 and replace the pick with the third character with probability 1/3. The same process is continued until the final character. In other words, let $C_n$ be the character that we pick after we have scanned $n$ characters and the $(n+1)$th character exists, the probability of keeping the pick is $\frac{n}{n+1}$ and the probability of switching to the $(n+1)$th character is $\frac1{n+1}$. Using simple induction, we can easily prove that each character has $1/m$ probability of being chosen if there are $m$ characters.

Search algorithm

A. Develop an algorithm to find both the minimum and the maximum of $n$ numbers using no more than $3n/2$ comparisons.

Solution: For an unsorted array of $n$ numbers, it takes $n-1$ comparisons to identify either the minimum or the maximum of the array. However, it takes at most $3n/2$ comparisons to identify both the minimum and the maximum. If we separate the elements to $n/2$ pairs, compare the elements in each pair and put the smaller one in group 178 $A$ and the larger one in group $B$. This step takes $n/2$ comparisons. Since the minimum of the whole array must be in group $A$ and the maximum must be in group $B$, we only need to find the minimum in $A$ and the maximum in $B$, either of which takes $n/2-1$ comparisons. So the total number of comparisons is at most $3n/2$.⁵

B. You are given an array of numbers. From the beginning of the array to some position, all elements are zero; after that position, all elements are nonzero. If you don’t know the size of the array, how do you find the position of the first nonzero element?

Solution: We can start with the 1st element; if it is zero, we check the 2nd element; if the 2nd element is zero, we check the 4th element… The process is repeated until the $i$th step when the $2^i$th element is nonzero. Then we check the $\frac{2^i+2^{i-1}}2$th element. If it is zero, the search range is limited to the elements between the $\frac{2^i+2^{i-1}}2$th element and the $2^i$th element; otherwise the search range is limited to the elements between the $2^{i-1}$th element and the $\frac{2^i+2^{i-1}}2$th element… Each time, we cut the range by half. This method is basically a binary search. If the first nonzero element is at position $n$, the algorithm complexity is $\Theta(\log n)$.

C. You have a square grid of numbers. The numbers in each row increase from left to right. The numbers in each column increase from top to bottom. Design an algorithm to find a given number from the grid. What is the complexity of your algorithm?

Solution: Let $A$ be an $n\times n$ matrix representing the grid of numbers and $x$ be the number we want to find in the grid. Begin the search with the last column from top to bottom: $A_{1,n},\cdots,A_{n,n}$. If the number is found, then stop the search. If $A_{n,n}\lt x$, $x$ is not in the grid and the search stops as well. If $A_{i,n}\lt x\lt A_{i+1,n}$, then we know that all the numbers in rows $1,\cdots,i$ are less than $x$ and are eliminated as well.⁶ Then we search the $(i+1)$th row from right to left. If the number is found in the $(i+1)$th row, the search stops. If $A_{1,i+1}>x$, $x$ is not in the grid since all the number in rows $i+1$ and above are larger than $x$. If $A_{i+1,j+1}>x>A_{i+1,j}$, we eliminate all the numbers in columns $j+1,\cdots,n$. Then we can search along column from $A_{i+1,j}$ towards $A_{n,j}$ until we find $x$ (or $x$ does not exist in 179 the grid) or a $k$ that makes $A_{k,j}\lt x\lt A_{k+1,j}$ and then we search left along the row $k+1$ from $A_{k+1,j}$ towards $A_{k+1,1}$… Using this algorithm, the search takes at most $2n$ steps. So its complexity is $O(n)$.

Fibonacci numbers

Consider the following C++ program for producing Fibonacci numbers:

int Fibonacci(int n)
{
    if (n <= 0)
        return 0;
    else if (n==1)
        return 1;
    else
        return Fibonacci(n-1)+Fibonacci(n-2);
}

If for some large $n$, it takes 100 seconds to compute Fibonacci($n$), how long will it take to compute Fibonacci($n+1$), to the nearest second? Is this algorithm efficient? How would you calculate Fibonacci numbers?

Solution: This C++ function uses a rather inefficient recursive method to calculate Fibonacci numbers. Fibonacci numbers are defined as the following recurrence:

$F_n$ has closed-formed solution $\displaystyle F_n=\frac{(1+\sqrt5)^n-(1-\sqrt5)^n}{2^n\sqrt5}$, which can be easily proven using induction. From the function, it is clear that

So the running time is a proportional to a sequence of Fibonacci numbers as well. For a large $n$, $(1-\sqrt5)^n\to0$, so $\displaystyle\frac{T(n+1)}{T(n)}\approx\frac{\sqrt5+1}{2}$. If it takes 100 seconds to compute Fibonacci($n$), the time to compute Fibonacci($n+1$) is $\displaystyle T(n+1)\approx\frac{\sqrt5+1}{2}T(n)\approx162$ seconds.⁷

The recursive algorithm has exponential complexity $\Theta(((\sqrt5+1)/2)^n)$, which is surely inefficient. The reason is that it fails to effectively use the information from Fibonacci numbers with smaller $n$ in the Fibonacci number sequence. If we compute $F_0,F_1,\cdots,F_n$ in sequence using the definition, the running time has complexity $\Theta(n)$.

An algorithm called recursive squaring can further reduce the complexity to $\Theta(\log n)$.

Since $\displaystyle\begin{bmatrix}F_{n+1}&F_n\\F_n&F_{n-1}\end{bmatrix}=\begin{bmatrix}1&1\\1&0\end{bmatrix}\times\begin{bmatrix}F_n&F_{n-1}\\F_{n-1}&F_{n-2}\end{bmatrix}$ and $\displaystyle\begin{bmatrix}F_2&F_1\\F_1&F_0\end{bmatrix}=\begin{bmatrix}1&1\\1&0\end{bmatrix}$, we can show that $\displaystyle\begin{bmatrix}F_{n+1}&F_n\\F_n&F_{n-1}\end{bmatrix}=\begin{bmatrix}1&1\\1&0\end{bmatrix}^n$ using induction. Let $\displaystyle A=\begin{bmatrix}1&1\\1&0\end{bmatrix}$, we can again apply the divide-and-conquer paradigm to calculate $A^n$:

The multiplication of two $2\times2$ matrices has complexity $\Theta(1)$. So $T(n)=T(n/2)+\Theta(1)$. Applying the master theorem, we have $T(n)=\Theta(\log n)$.

Maximum contiguous subarray

Suppose you have a one-dimensional array $A$ with length $n$ that contains both positive and negative numbers. Design an algorithm to find the maximum sum of any contiguous subarray $A[i,j]$ of $A$: $\displaystyle V(i,j)=\sum_{x=i}^jA[x],\ 1\le i\le j\le n$.

Solution: Almost all trading systems need such an algorithm to calculate maximum run-up or maximum drawdown of either real trading books or simulated strategies. Therefore this is a favorite algorithm question of interviewers, especially interviewers at hedge funds and trading desks.

The most apparent algorithm is an $O(n^2)$ algorithm that sequentially calculates the $V(i,j)$’s from scratch using the following equations:

$V(i,i)=A[i]$ when $j=i$ and $\displaystyle V(i,j)=\sum_{x=i}^jA[x]=V(i,j-1)+A[j]$ when $j>i$.

As the $V(i,j)$’s are calculated, we also keep track of the maximum of $V(i,j)$ as well as the corresponding subarray indices $i$ and $j$.

A more efficient approach uses the divide-and-conquer paradigm. Let’s define $\displaystyle T(i)=\sum_{x=1}^iA[x]$ and $T(0)=0$, then $V(i,j)=T(j)-T(i-1),\ \forall1\le i\le j\le n$. Clearly for any fixed $j$, when $T(i-1)$ is minimized, $V(i,j)$ is maximized. So the maximum subarray ending at $j$ is $V_{\max}=T(j)-T_{\min}$ where $T_{\min}=\min(T(1),\cdots,T(j-1))$. If we keep track of and update $V_{\max}$ and $T_{\min}$ as $j$ increases, we can develop the following $O(n)$ algorithm:

T = A[1]; Vmax = A[1]; Tmin = min(0,T)
For j = 2 to n
{ T = T + A[j];
    If T-Tmin > Vmax then Vmax = T-Tmin;
    If T < Tmin, then Tmin = T;
}
Return Vmax;

The following is a corresponding C++ function that returns $V_{\max}$ and indices $i$ and $j$ given an array and its length:

double maxSubarray(double A[], int len, int &i, int &j)
{
    double T=A[0], Vmax=A[0];
    double Tmin = min(0.0, T);
    for(int k=1; k<len; ++k)
    {
        T+=A[k];
        if (T-Tmin > Vmax) {Vmax=T-Tmin; j=k;}
        if (T<Tmin) {Tmin = T; i = (k+1<j)? (k+1):j;}
    }
    return Vmax;
}

Applying it to the following array $A$,

double A[]={1.0,2.0,-5.0,4.0,-3.0, 2.0, 6.0, -5.0, -1.0};
int i = 0, j =0;

double Vmax = maxSubarray(A, sizeof(a)/sizeof(A[1]), i, j);

will give $V_{\max}=9$, $i=3$ and $j=6$. So the subarray is [4.0, −3.0, 2.0, 6.0].