Generating 1/f Noise Sequences as Constraint Satisfaction: The Voss Constraint

1 1/f noise in nature and human artefacts

Many natural phenomena have been shown to exhibit socalled 1/f fluctuations, such as fluctuations of river flooding or luminosity of stars [Mandelbrot, 1982]. Such fluctuations have also been observed in human artefacts, such as fluctuations in rhythm [Hennig et al., 2011], or music [Voss and Clarke, 1975]. 1/f fluctuations have also been shown to be ubiquitous in human cognition (see, e.g., [Farrell et al., 2006a; 2006b]).

1/f noise, also called pink noise, is defined by a power law relation between the frequency f of a signal and its power spectral density S(f) :

The ubiquity of 1/f fluctuations has therefore been used to enhance artificial artefacts, to make them look or sound more natural or more human. For instance, [Hennig et al., 2011] introduce 1/f fluctuations in computer-generated musical rhythms to make them sound less artificial, and show experimentally that these subtle differences are indeed perceived and preferred by humans. Similarly, 1/f noise is often used to make images or textures more natural, see e.g., Perlin noise [Perlin, 1985]. Many algorithms have been designed to generate pink noise series (see, e.g., [Kasdin, 1995]), but these algorithms are always designed as self-contained black boxes. As a consequence, it is difficult to generate 1/f series that satisfy simultaneously additional properties.

The striking result obtained by [Voss and Clarke, 1978] is that 1/f pitch sequences sound more natural than sequences obtained with other distributions (random, Brownian). However, most music examples of 1/f sequences exhibited in the literature do not satisfy other basic properties of music such as meter, harmony, or higher-level structure. For instance, the melodies generated in [Voss and Clarke, 1978] do not contain bars, i.e., are not metrically correct: 1/f noise is a global property that is not easily combined with other properties. In that case, generating sequences that fulfill metrical properties is in itself a difficult combinatorial problem [Roy and Pachet, 2013]. In the context of sequence generation, it is therefore natural to look for an 1/f global constraint, which can be added to other arbitrary constraints, in the framework of constraint satisfaction. The goal of this paper is to propose a simple way to do so in the context of hard constraints, by taking inspiration from an algorithm proposed by physicist Richard Voss and popularized by Martin Gardner in [Gardner, 1978]. Like Voss' algorithm, the scheme we propose does not solve the problem in its full generality, but provides a simple and elegant solution to the issue of 1/f sequence generation. Before describing this algorithm and its formulation as a constraint satisfaction problem, we review the state-of-the-art in global constraints dealing with distributions.

Algorithm 1: Voss' algorithm as described by Gardner.

Input: N a sequence length, nbDice a number of dice, maxDice the maximum dice value Output: A sequence of N values in [nbDice, nbDice.maxDice] with a 1/f spectrum result ← new integer array of length N diceValues ← new integer array of length nbDice for i = 0, . . . , N - 1 do for d = 1, . . . , nbDice do if (dth bit of i) ̸= (dth bit of i - 1) then // toss die number d diceValues (d) ← Random integer in [1,maxDice] // sum all dice result (i) ← ∑nbDice k=1 diceValues(k) return result

2 Related work

The spread constraint [Pesant and R´egin, 2005] ensures that the values in a sequence have a given mean and standard deviation. Similarly, the deviation constraint [Schaus et al., 2007] ensures that values of a set of variables have a given standard deviation. The balance constraint [Beldiceanu et al., 2007] bounds the difference in the number of occurrences of values assigned to variables. No constraint, to our knowledge, deals with more complex distribution of values, in particular taking into account the spectrum of the sequence (seen as a time series). Indeed, the spectrum of a time series is obtained through a complex operation (the Fourier transform), and, to our knowledge, no Spectrum constraint has yet been proposed. Recently, the idea of encoding statistical systems as constraint optimization problems has emerged, notably for Markov chains [Morin and Quimper, 2014; Pachet and Roy, 2011] and neural networks [Lombardi and Gualandi, 2013; Bartolini et al., 2011]. This work can be seen as yet another bridge between stochastic algorithms and finite-domain, discrete constraint satisfaction.

3 Voss' algorithm

In our case we do not need to control the spectrum of a sequence precisely but only to enforce a specific, global property of the spectrum. The seminal paper [Gardner, 1978] described a simple algorithm to generate 1/f sequences, i.e., sequences of integers having a spectral density in 1/f. This algorithm, invented by Voss but never published to our knowledge, consists in using a number of dice and summing their values, exploiting the fact that 1/f series can be modeled as sums of random variables with specific correlations. At each step, only certain dice are tossed: those who correspond to a bit change when writing a sequence of integers in base 2. As a consequence, it can be seen easily that in the resulting sequence the variance is inversely proportional to the frequency, due to the fact that some dice stay unchanged over longer periods than others. The algorithm is described in Algorithm 1 and the dice tossing scheme is illustrated in Table 1. Three 6-sided dice are considered, represented by 3 colors (Blue, Green, Red) to generate a sequence of 8 values. Each line represents an integer, from 0 to 7 written in base 2. The numbers of the sequence are obtained by summing up the values of the three dice at each line. The trick is that only dice corresponding to a bit change with respect to the previous line are tossed.

Voss' algorithm generates 1/f series, or at least a very good approximation thereof. We can show it by estimating the power spectral density using S(f) = 1/T |X(f, T)|2 , where X(f, T) = ∫ T 0 x(t)e-iftdt is the Fourier transform of the sequence. We then perform a simple linear regression in the log-log spectrum in order to obtain the slope -α. Such a series is illustrated in Figure 1, and its log-log spectrum in Figure 2, which is well approximated by a straight line, at least for values larger than a threshold (-2.5).

4 The Voss constraint

Voss' algorithm is a statistical one, and cannot be, as such, translated easily into a hard constraint satisfaction problem. However, the decomposition of 1/f series as sums of random variables can lead to a simple representation as a set of constraints, in combination with a random heuristics for variable value choice. The idea is to consider both sequence elements and dice as constrained variables and relate them through a set of sum constraints [Trick, 2003].

The Voss constraint can be posted on sequences of size N = 2K , whose elements have a domain of the form D = [K,K.R], with R a positive integer. The Voss constraint ensures that the values of the variables could have been obtained by tossing dice according to Voss' algorithm presented above, with parameters N, K for nbDice and R for maxDice.

We first present a naive implementation of the Voss constraint, and then a more efficient one.

4.1 A naive CSP

The Voss constraint can be implemented naively by introducing a die variable representing each die toss. Die variable Di,j represents the j th toss of die i (first index is die number and second one is toss index). We then relate the sequence variables to sums of these die variables, according to Voss' algorithm.

More precisely, we introduce the Voss tree that represents the structure of the die variables as follows.

Definition 1 ((K,R)-Voss tree). For a given integer K and R, the (K,R)-Voss tree is a binary tree of variables Di,j, with i = 1, . . . ,K and j = 1, . . . , 2i . Each Di,j has two children Di+1,2j-1 and Di+1,2j. The domain of each Di,j is {1, . . . , R}. K is called the order of the tree, R is its range.

The total number of variables in a (K,R)-Voss tree is

which is less than 2×2K .

For each sequence variable Xi we introduce a sum constraint equating Xi with the sum of the die variables obtained by walking up the Voss tree from the leaf Dk,i:

Figure 3 shows the Voss tree for a sequence of length 8, thus K = 3. In that case, the number of die variables is 14, and we post the 8 constraints shown in Figure 4.

For a sequence of length N = 2K the Voss tree contains 2N - 2 extra die variables (Equation 1).

X1 = D1,1 +D2,1 +D3,1,

X2 = D1,1 +D2,1 +D3,2,

. . .

X8 = D1,2 +D2,4 +D3,8.

Figure 4: The n-ary sum equations defining the Voss constraint.

Definition 2 ((K,R)-Voss constraint). Let X1, . . . ,X2K be a sequence of variables, each with finite domain D = [K,K.R]. The Voss constraint is the conjunction of all sum constraints of the (K, R)-Voss tree defined by Equation 2.

Filtering the Voss constraint can be achieved by filtering each sum constraint individually. However, domainconsistency of sum constraints is pseudo-polynomial [Trick, 2003]. Additionally, domain-consistency for the whole Voss constraint is not guaranteed by domain consistency of individual sum constraints.

4.2 An efficient filtering for the Voss constraint

A much more efficient filtering method for the Voss constraint is obtained by replacing each n-ary sum constraint by a set of ternary constraints representing partial sums. More precisely, for each Di,j we introduce a variable Si,j representing the partial sum from Di,j to the top of the Voss tree. In order to avoid computing n-ary sums, Si,j are defined recursively with ternary sum constraints as follows:

The sequence variables are given by Xi = SK,i.

The number of sum variables (Si,j) is equal to the number of die variables so the total number of extra variables (Di,j and Si,j) is 2K+2 - 4, also linear in the size of the sequence. The number of ternary sum constraints is 2K+1 - 2. The ternary constraints corresponding to our running example (length 8, 3 dice) are illustrated in Figure 5 and Figure 6.

By construction, the resulting CSP is Berge-acyclic, i.e., the network is a tree and the constraints overlap on, at most, one variable. Consequently, domain consistency for the whole constraint is achieved by propagating each ternary sum

Figure 6: The ternary sum equations defining the Voss constraint.

constraint individually [Beeri et al., 1983; J´egou, 1993]. If sequence variables have a domain size d, each additional variable Si,j has a domain bounded by d. The complexity of filtering one such ternary constraint is in O(d2 ). Therefore, for a sequence of length N, the complexity of achieving domainconsistency for the Voss constraint is in O(N.d2 ).

4.3 Generalization

The constraint can be easily generalized in two ways:

• Arbitrary sequence length. The running example focuses on the case where the length is a power of 2. The tree structure proposed can easily accommodate arbitrary sequence length. The resulting tree is incomplete but the incidence on the spectrum is negligible. For a given sequence length l, the number of dice to consider is nbDice = ⌈log 2(l)⌉.
• Arbitrary domains sizes. The range of dice maxDice or R is obviously related to the domain size, since values are by construction in [N,N ∗ R]. For a contiguous domain [min,max], R is to be chosen as [⌊min/nbDice⌋, ⌈max/nbDice⌉].

5 AMelody Generation Example

We revisit the classical melody generation example introduced by Gardner by generating again melodies. Our goal is not to show the effect of 1/f sequences on the perception of melodies, which has already been done convincingly in [Voss and Clarke, 1978]. Our goal is to show that 1) the Voss constraint does generate sequences with a 1/f spectral density, like the original Voss algorithm, and 2) that this constraint, when added to a CSP with another constraint, here a GCC, still biases the spectrum of the resulting sequence to be in 1/f .

5.1 Implementation

We implemented the Voss constraint in BackJava, an inhouse Java finite-domain constraint solver similar in nature to Choco [choco Team, 2010], as well as the MusES musical object library [Pachet, 1994]. All the experiments ran on a machine with a Core i7, 2.3 GHz CPU, with 16GB RAM, and running an Oracle Java 7 VM under Windows 8.

5.2 Generating LongMelodies

We consider a melody generation problem defined by three sets of constraints. We generate sequences of integers (representing the pitch of musical notes) of length 512, so that the spectrum can be computed without side effects. The range of pitches (i.e., the domain D of sequence variables ) is [0, 16], i.e., covers roughly 2 octaves in a diatonic setting. The number of dice is 9 (= log2(512)), and their range is [0, 1, 2]. We consider three cases:

1. A Voss constraint only holding on the whole sequence,
2. A global cardinality constraint (GCC) [R´egin, 1996] only. The GCC ensures that there is a balanced distribution of values in the resulting sequence. We post the GCC as follows: We consider Nopt the optimum distribution, that would ensure that all values are equally represented: Nopt = N/ |D|. We enforce a GCC forcing each value of D to occur exactly this optimum value,
3. A Voss constraint and the GCC constraint.

For each case we compute the log-log spectrum and estimate the slope of curve as previously. The spectrum of case #1 (Voss constraint only, see Figure 11) is shown in Figure 7. It can be seen clearly that the spectrum is in 1/f. For case #2 (GCC only, Figure 12) the spectrum (Figure 8) is clearly not in 1/f. Case #3 (same GCC plus the Voss constraint, see Figure 13) shows that the spectrum (Figure 9) is again in 1/f thanks to the Voss constraint. The distribution of the slopes of the various spectra for 100 solutions is shown in Figure 10 and confirms that the presence of the Voss constraint does bias the spectrum toward an 1/f distribution. Running times and number of backtracks, averaged for one solution, are reported for each setup on Table 2. It can be seen that adding a Voss constraint increases the cost of finding solutions in a reasonable way.

5.3 Shorter Melodies

The impact of the Voss constraint on shorter melodies can be assessed through examples rather than from the analysis of the spectrum. Here, melodies are generated from a Markov model estimated from two well-known songs of Michel Legrand ("You must believe in Spring" and "How do you keep the music playing?").

In the first case (Figure 14) we use the meter constraint [Roy and Pachet, 2013] to enforce metrical properties of the generated sequence (notes do not cross bars), and a total duration of 12 bars (in 4/4). Without the Voss constraint, meter has a tendency to generate solutions which include long cycles (see, e.g., the series of 5 consecutive C at the end of the melody in Figure 14). Also, meter (which uses a Markov model of order 1) may generate series of large intervals in a row, which sounds awkward. These two cases seem to be less frequent when using the Voss constraint, because of the imposition of the specific correlation structure (Figure 15).

6 Discussion

We have introduced the Voss constraint, which ensures that the elements of a sequence are generated from the dice tossing scheme invented by Voss as described by Gardner. We showed that the constraint can be implemented as a tree of ternary sum constraints, yielding an efficient domainconsistency procedure, so its overhead is arguably small. The Voss constraint can be added to any sequence generation problem to introduce a pressure for the sequence to exhibit a 1/f spectral density in average, and therefore look or sound more natural.

There are several limitations to our proposal. First, Voss' algorithm imposes a specific structure of correlation, corresponding to the ordering imposed by an increasing integer sequence written in base 2. However, any other ordering could be used as well. This can be addressed by a generalization of the algorithm in which the number of dice to toss is chosen at random at each step, which would also result in a 1/f distribution. The probability to toss d dice is proportional to 1/2d to match on average the frequency of dice re-tossing performed by Voss's original algorithm. To do so, however, requires the use of dynamic constraints (sums constraints would be posted depending on the value of this variable) or reified constraints, which would increase the filtering cost. Another approach is to draw a correlation structure at random before the resolution, but this does not ensure that the chosen structure is compatible with the other constraints.

Another extension of the algorithm concerns the value of α. Voss' algorithm produces exactly 1/f series (α = 1). It could be interesting to control the value of α, typically within [0,2]. This could be done, again, by introducing a variable representing the number of dice to toss at each step, and biasing its probability distribution through a specific value ordering heuristics. Another limitation concerns the coding. The Voss constraint requires a coding from the domain (e.g., musical notes) to integers (in our example the pitch of a note). This coding is arbitrary and the user may want to experiment with various codings (e.g., integers representing the degree of the note in a scale) and there is no best coding a priori. Finally, one could ensure that solutions look more random by enforcing additional constraints on dice. For instance one could enforce a uniformly distributed set of values, to simulate the effect of the law of big numbers. This could be done by posting cardinality or nValue [Bessi`ere et al., 2006] constraints on each die (i.e., for a die i, all the Di,j).

There are obviously limits to representing stochastic properties with hard constraints. Another approach, in progress, is to model 1/f as a statistical constraint in the line of Walsh [2002] or Rossi et al. [2014] using the same die structure. However, we think that the work presented here contributes to constraint-based sequence generation by providing insights on this elusive but ubiquitous 1/f phenomenon, seen from a generative perspective.

Acknowledgements

This research is conducted within the Flow Machines project which received funding from the European Research Council under the European Union's Seventh Framework Programme (FP/2007-2013) / ERC Grant Agreement n. 291156. We wish to thank Mirko Degli Esposti for insightful discussions about 1/f distributions.