Prefix sum array and difference array

From PEGWiki
Revision as of 08:16, 19 December 2011 by Brian (Talk | contribs) (rewrite --- should be easier to read now)

Jump to: navigation, search

Given an array of numbers, we can construct a new array by replacing each element by the difference between itself and the previous element, except for the first element, which we simply ignore. This is called the difference array, because it contains the first differences of the original array. We will denote the difference array of array A by D(A). For example, the difference array of A = [9, 2, 6, 3, 1, 5, 0, 7] is D(A) = [2-9, 6-2, 3-6, 1-3, 5-1, 0-5, 7-0], or [-7, 4, -3, -2, 4, -5, 7].

We see that the difference array can be computed in linear time from the original array, and is shorter than the original array by one element. Here are implementations in C and Haskell. (Note that the Haskell implementation actually takes a list, not an array, and will work fine for infinite lists.)

// D must have enough space for n-1 ints
void difference_array(int* A, int n, int* D)
{
    for (int i = 0; i < n-1; i++)
        D[i] = A[i+1] - A[i];
}

d :: [Int] -> [Int] -- NB: The general type signature should be (Num a) => [a] -> [a]
d a = zipWith (-) (tail a) a

The prefix sum array is the opposite of the difference array. Given an array of numbers A and an arbitrary constant c, we first append c onto the front of the array, and then replace each element with the sum of itself and all the elements preceding it. For example, if we start with A = [9, 2, 6, 3, 1, 5, 0, 7], and choose to append the arbitrary value -8 to the front, we obtain P(-8, A) = [-8, -8+9, -8+9+2, -8+9+2+6, ..., -8+9+2+6+3+1+5+0+7], or [-8, 1, 3, 9, 12, 13, 18, 18, 25]. Computing the prefix sum array can be done in linear time as well, and the prefix sum array is longer than the original array by one element:

// P must have enough space for n+1 ints
void prefix_sum_array(int c, int* A, int n, int* P)
{
    P[0] = c;
    for (int i = 0; i < n; i++)
        P[i+1] = P[i] + A[i];
}

p :: Int -> [Int] -> [Int] -- NB: The general type signature should be (Num a) => a -> [a] -> [a]
p = scanl (+)              -- NB: This is the same as p c a = scanl (+) c a

Note that every array has an infinite number of possible prefix sum arrays, since we can choose whatever value we want for c. For convenience, we usually choose c = 0. However, changing the value of c has only the effect of shifting all the elements of P(c,A) by a constant. For example, P(15, A) = [15, 24, 26, 32, 35, 36, 41, 41, 48]. However, each element of P(15, A) is exactly 23 more than the corresponding element from P(-8, A).

The functions D and P carry out reverse processes. Given an nonempty zero-indexed array A:

  1. D(P(c, A)) = A for any c. For example, taking the difference array of P(-8, A) = [-8, 1, 3, 9, 12, 13, 18, 18, 25] gives [9, 2, 6, 3, 1, 5, 0, 7], that is, it restores the original array A.
  2. P(A_0, D(A)) = A. Thus, taking D(A) = [-7, 4, -3, -2, 4, -5, 7] and A_0 = 9 (initial element of A), we have P(A_0, D(A)) = [9, 2, 6, 3, 1, 5, 0, 7], again restoring the original array A.

Analogy with calculus

These two processes—computing the difference array, and computing a prefix sum array—are the discrete equivalents of differentiation and integration in calculus, which operate on continuous domains. An entry in an array is like the value of a function at a particular point.

  • Reverse processes:
  • D(P(c, A)) = A for any c. Likewise \frac{d}{dx} \int_c^x f(t)\, dt = f(x) for any c.
  • P(A_0, D(A)) = A. Likewise f(a) + \int_a^x \frac{df}{dt}\, dt = f(x).
  • Uniqueness:
  • A differentiable function f(x) can only have one derivative, \frac{df}{dx}. An array A can only have one difference array, D(A).
  • A continuous function f(x) has an infinite number of antiderivatives, F_c(x) = \int_c^x f(t)\, dt, where c can be any number in its domain, but they differ only by a constant (their graphs are vertical translations of each other). An array A has an infinite number of prefix arrays P(c,A), but they differ only by a constant (at each entry).
  • Given some function f:[a,b]\to\mathbb{R}, and the fact that F, an antiderivative of f, satisfies F(a) = y_0, we can uniquely reconstruct F. That is, even though f has an infinite number of antiderivatives, we can pin it down to one once we are given the value the antiderivative is supposed to attain on the left edge of f's domain. Likewise, given some array A and the fact that P, a prefix sum array of A, satisfies P_0 = c, we can uniquely reconstruct P.
  • Effect on length:
  • D(A) is shorter than A by one element. Differentiating f:[a,b] \to \mathbb{R} gives a function f':(a,b) \to \mathbb{R} (shortens the closed interval to an open interval).
  • P(c,A) is longer than A by one element. Integrating f:(a,b) \to \mathbb{R} gives a function F:[a,b] \to \mathbb{R} (lengthens the open interval to a closed interval).

Because of these similarities, we will speak simply of differentiating and integrating arrays. An array can be differentiated multiple times, but eventually it will shrink to length 0. An array can be integrated any number of times.

Use of prefix sum array

The Fundamental Theorem of Calculus also has an analogue, which is why the prefix sum array is so useful. To compute an integral \int_a^b f(t)\, dt, which is like a continuous kind of sum of an infinite number of function values f(a), f(a+\epsilon), f(a+2\epsilon), ..., f(b), we take any antiderivative F, and compute F(b) - F(a). Likewise, to compute the sum of values A_i, A_{i+1}, A_{i+2}, ..., A_{j-1}, we will take any prefix array P(c,A) and compute P_j - P_i. Notice that just as we can use any antiderivative F because the constant cancels out, we can use any prefix sum array because the initial value cancels out. (Note our use of the left half-open interval.)

Proof: P_j = c + \sum_{k=0}^{j-1} A_k and P_i = c + \sum_{k=0}^{i-1} A_k. Subtracting gives P_j - P_i = \sum_{k=0}^{j-1} A_k - \sum_{k=0}^{i-1} A_k = \sum_{k=i}^{j-1} A_k as desired. _{\blacksquare}

This is best illustrated via example. Let A = [9,2,6,3,1,5,0,7] as before. Take P(0,A) = [0, 9, 11, 17, 20, 21, 26, 26, 33]. Then, suppose we want A_2 + A_3 + A_4 + A_5 = 6 + 3 + 1 + 5 = 15. We can compute this by taking P_6 - P_2 = 26 - 11 = 15. This is because P_6 - P_2 = (0 + A_0 + A_1 + A_2 + A_3 + A_4 + A_5) - (0 + A_0 + A_1) = A_2 + A_3 + A_4 + A_5.

Example: Counting Subsequences (SPOJ)

Computing the prefix sum array is rarely the most difficult part of a problem. Instead, the prefix sum array is kept on hand because the algorithm to solve the problem makes frequent reference to range sums.

We will consider the problem Counting Subsequences from IPSC 2006. Here we are given an array of integers S and asked to find the number of contiguous subsequences of the array that sum to 47.

To solve this, we will first transform array S into its prefix sum array P(0,S). Notice that the sum of each contiguous subsequence S_i + S_{i+1} + S_{i+2} + ... + S_{j-1} corresponds to the difference of two elements of P, that is, P_j - P_i. So what we want to find is the number of pairs (i,j) with P_j - P_i = 47 and i < j. (Note that if i > j, we will instead get a subsequence with sum -47.)

However, this is quite easy to do. We sweep through P from left to right, keeping a map of all elements of P we've seen so far, along with their frequencies; and for each element P_j we count the number of times P_j - 48 has appeared so far, by looking up that value in our map; this tells us how many contiguous subsequences ending at S_{j-1} have sum 47. And finally, adding the number of contiguous subsequences with sum 47 ending at each entry of S gives the total number of such subsequences in the array. Total time taken is O(N), if we use a hash table implementation of the map.