Multiple scaling behavior and nonlinear traits in music scores

We present a statistical analysis of music scores from different composers using detrended fluctuation analysis. We find different fluctuation profiles that correspond to distinct auto-correlation structures of the musical pieces. Further, we reveal evidence for the presence of nonlinear auto-correlations by estimating the detrended fluctuation analysis of the magnitude series, a result validated by a corresponding study of appropriate surrogate data. The amount and the character of nonlinear correlations vary from one composer to another. Finally, we performed a simple experiment in order to evaluate the pleasantness of the musical surrogate pieces in comparison with the original music and find that nonlinear correlations could play an important role in the aesthetic perception of a musical piece.


Introduction
Music is a complex construct that involves cultural factors, acoustic features, interpretation techniques and audience perception. Its ample range of properties make it a fascinating study subject from many different viewpoints, in addition to its relevance in our everyday life. The study of music from the perspective of statistical physics has been of great interest during the last decades. One of the pioneering works in this context, using techniques from statistical physics, was published by Voss & Clarke [1]. They estimated the spectral density of intensity fluctuations and found a power law behaviour close to 1/f noise. This result inspired many researchers to study the statistical properties of music, ranging from the identification of temporal patterns or power laws, to the development of algorithms for music composition [2][3][4][5][6][7][8][9][10][11]. Among these, are the identification of mayor and minor tonalities in Bach's welltempered clavier by means of a statistical parametrization [6], 2017 The Authors. Published by the Royal Society under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/4.0/, which permits unrestricted use, provided the original author and source are credited.

Construction of the time series
We consider music scores as a sequence of integer numbers, each of them labelling a different note. The numbers were extracted from midi files obtained from different web databases [17,18]. We processed the midi files using midicsv [19], a free software that converts midi into csv (comma-separated values) files. Taking the note with the smallest duration as the time unit, a time series can be generated as shown in figure 1. An improved understanding of such representation can be attained from figures 2 and 3. In figure 3 the first eight measures of the piece shown in figure 2 are translated in units of the shortest note (in the figure we subdivided each measure into eighth notes, which represent the shortest duration of the original piece and serves as the time unit in this case). The time series is multivariate and the number of variables depends on the number of instruments or voices of the piece.

Detrended fluctuation analysis
The detrended fluctuation analysis (DFA) method, developed by Peng et al. [20] (see a detailed and pedagogic presentation in [21]), was introduced in order to avoid the detection of spurious correlations generated by trends in time series. It has been used in other works of music analysis providing promising results [4,9,10]. As it is unclear whether a given musical piece is stationary or not, DFA seems to be a suitable method for a proper measurement of correlations.    The procedure is repeated by varying s such that the fluctuation function is obtained in terms of the segment length, which represents the timescale where correlations might be present. When autocorrelations scale like a power law, the RMS fluctuation function F(s) behaves as F(s) ∼ s α , where α is the Hurst exponent. A value of α > 0.5 indicates the presence of persistent correlations, e.g. α = 1 is the case for 1/f noise. On the other hand, a value of 0 < α < 0.5 corresponds to anticorrelations and α = 0.5 to white noise [20]. Music scores, in general, can be represented as multivariate datasets where to each voice (or instrument, respectively) a data channel is assigned. For multivariate series with n dimension we compute the n-DFA [10,22] with the RMS function: where z(j) is a vector whose components contain the value of the pitch of each voice from the original score at time point j and z m the same for the polynomial fits to each voice. If one observes a power law for the fluctuation function, one can relate the Hurst exponent α to the exponent of the power spectrum

Magnitude and sign detrended fluctuation analysis
Ashkenazy et al. [24,25] developed a variation of the DFA method capable to detect nonlinear autocorrelations within empirical recordings. This method can be summarized by the following recipe: (i) for a given time series x(i) the increment series is defined as (ii) the increment series is decomposed into a magnitude series and sign series: x(i) = sgn( x(i))| x(i)|, their respective means are subtracted to avoid artificial trends; (iii) because of the limitations of the DFA method for estimating α < 0.5 (anti-correlated series), the magnitude and sign series are integrated first to make sure they are positively correlated [25]. (iv) The DFA method is implemented on the integrated magnitude and sign series. (v) To obtain the respective scaling exponents, the function F(s)/s is estimated; the 1/s factor is to compensate for the integration made before. If the data obey a scaling law, the fluctuation function should behave as F(s)/s ∼ s α−1 . It has been shown that the magnitude series carries information regarding nonlinear properties of the original time series [24]. All DFA estimations presented in this study are multivariate and are performed with a polynomial of degree 2 (n-DFA-2). For simplicity, we will refer to them merely as DFA. We also tried higherorder polynomials m = 3, 4, obtaining quantitatively similar results. We further denote the second-order multivariate magnitude DFA (n-MDFA-2) just as MDFA.

Surrogate data
To validate the significance of the results obtained by the magnitude DFA, it is necessary to compare them with those obtained for appropriately generated surrogate data, which represent the nullhypothesis of zero nonlinear autocorrelations. Hence, while generating surrogates from the original data, one should destroy all nonlinear features which the original data may contain, while conserving the same linear autocorrelations. The complete information about the linear correlation structure is imprinted in the power spectral density (the distribution of the square of the amplitudes of the complex Fourier coefficients). Nonlinear correlations, on the other hand, are inherent in the distribution of the Fourier phases. In this study, surrogate data are generated in an iterative fashion where the amplitude distribution as well as the power spectrum are adjusted to those of the original data, while Fourier phases are replaced by random numbers, uniformly distributed between zero and 2π [26]. These time series share the same linear univariate properties as the original recordings, but lack their nonlinear correlations. All surrogate data used in this study were generated with the freely available TISEAN package. [27] 3. Results

Linear correlations
We applied both the DFA as well as the MDFA method to 304 music scores of different composers listed in table 1.
We find that not all of the DFA functions have a single power law scaling; some of them show a crossover from one scaling region to another, while others do not have any power law behaviour on any scale. For this reason, it is not possible to estimate a Hurst exponent for every music score. However, we find that the function log(F(s)) still provides interesting information about the autocorrelation structure of the pieces. We could identify five qualitatively different profiles in the fluctuation function, which rsos.royalsocietypublishing.org R. Soc. open sci.   characterize the autocorrelation structure of the musical piece. The profiles are shown in figure 4 and table 2.
The first profile, exemplified by a piece from Palestrina in figure 4 indicates strong short-range autocorrelations, which are getting weaker at longer timescales. In these pieces, the memory of certain motifs of the musical structure of small durations gets lost on longer timescales such that self-similarity is diluted and irregularity increases. Such type of crossovers have already been identified in [4].
Counterintuitively, the second profile shows opposite characteristics: correlations at large timescales are stronger than those at short timescales (see the example of Shostakovich in figure 4). These large time correlations, which are similar in magnitude to those of the first profile, are indicative of the recurrence of long patterns. However, in this case short sequences of notes are more irregular.
A third profile consists of a constant scaling behaviour over the whole range of box sizes s. An example is provided by a piece of Mozart in figure 4. The two remaining profiles correspond to different changes of curvature, illustrated by Beethoven and Bach in figure 4, showing a more irregular behaviour, where stronger and less strong autocorrelated timescales alternate. Table 2 summarizes the number of opuses of a given composer, which could be assigned to one of the profiles of the fluctuation function. Although table 2 is not necessarily conclusive, it provides some trends, which allow us to roughly distinguish composers. For instance, profile 1 (stronger shortrange than long-range correlations) and profile 3 (a clear scaling over the whole range) are the preferred correlation profiles for all composers except Dvorak, for whom we identified an important percentage of pieces with fluctuation profile 2 in comparison to profile 1. In particular, Palestrina shows a clear  Table 2. Classification by fluctuation profiles of the pieces analysed in this study. The first row shows the five different profiles we could identify. Numbers indicate the amount of scores of each composer with a given profile. The classification has been done by careful eye inspection of the scaling plots. preference for profile 1. On the other hand, pieces by Shostakovich, Dvorak and to some extent also those by Beethoven, are, relative to the composers selected in this study, most frequently assigned to profile 2. There, stronger autocorrelations act on longer timescales. Furthermore, while Beethoven inclines to profile 3, Mozart and Haydn show a preference for profile 1. Finally, Bach has highest scoring at profile 4 and Haydn for profile 5.
We find it interesting that the number of pieces corresponding to the second profile apparently increases the later the musical period of the composer. While the percentage of such a profile detected in this study is about 3 to 4% for Bach, Haydn and Mozart, this fraction increases to 11% for Beethoven and reaches values of approximately 36% and 30% for Dvorak and Shostakovich, respectively. Furthermore, it seems that Shostakovich shows the highest variability in terms of different autocorrelation structures. Weaker correlations at short distances might be a consequence of less restrictive composition rules at neighbouring notes, because an increase in the absolute values of the autocorrelation measures an increase in deviations from statistically independent fluctuations. The preference for the second profile clearly correlates with the period of time, and might be a direct cause of how the composition rules have changed in time.
A selection of DFA results for different composers, where a clear scaling behaviour over the whole range could be identified, is shown in figure 5. Estimates for the Hurst exponents are included in each graph. It is interesting that from Palestrina and Bach to the classical composers the exponent decreases, but then it increases drastically in Shostakovich. This behaviour is confirmed by the results presented in figure 6, which sums up the results for all cases with unique scaling. The largest Hurst exponents are obtained from musical pieces by Palestrina, whose cumulative distribution is well separated from the rest. This distribution is followed by that of Bach. The corresponding distribution is again somewhat separated from the pieces of Mozart, Beethoven and Haydn. On an average, the lowest values are found for the pieces composed by Dvorak. Then, the values for the pieces by Shostakovich are similar to those obtained for Bach. Hence, figure 6 is consistent with the observations of figure 5.
The higher linear autocorrelations in Palestrina and Bach can probably be attributed to the musical form and composition rules: stronger correlations are possibly related to the contrapuntal form from the Renaissance and Baroque periods, where theme repetitions and variations are present in every voice of the piece, a feature which contributes directly to the amount of linear autocorrelations.
The existence of a crossover in the DFA function of pitch sequences has been reported previously in [4] where palindromic compositions from Mozart are shown to have a separation of scaling properties at different ranges with α 1 > α 2 . The high recurrence of short temporal patterns is not extended to structures of longer duration, so correlations diminish after the crossover to large timescales. Some examples of a crossover in the fluctuation function are shown in figure 7 with α 1 > α 2 and α 2 > α 1 .
By simple eye revision, the presence of two different timescales can be appreciated in the patterns of the original time series shown in figures 8 and 9. Within the range of 1-20 notes in Bach's sonata we can observe many quite similar patterns of pitch fluctuations; this similarity weakens on larger timescales. This explains the larger slope log(F(s)) at the first timescale that decreases from α 1 = 1.4 to α 2 = 0.9 for distance s > 20, which is indicative of weaker correlations. The opposite happens in the case of the fugue from Shostakovich (figure 9), where the correlations on shorter timescales are considerably weak; they change from α 1 = 0.6 to α 2 = 0.9 at s ≈ 32. At short timescales the pitch fluctuations are more irregular and close to white noise (α = 0.5), hence melodic patterns are more difficult to predict. However, on longer timescales motives of note sequences are (almost) repeated, which explains the presence of stronger long-range correlations. Figure 10 shows a plot α 1 versus α 2 taking into account all musical pieces of the different composers which show a crossover in their fluctuation function. The centre of each ellipse is defined by the mean of each (α 1 , α 2 )-distribution, while the axes reflect their respective standard deviations. Note that the interval of variation of α 1 is considerably larger than the one for α 2 ; this is an indication that changes in long-range correlations are more restricted than in short-range ones. Thus, the overall structure of the musical forms is more preserved among the various composers than are the motifs.    in table 2, where Dvorak and Shostakovich showed a more frequent preference for the second autocorrelation profile with α 2 > α 1 . In general, we found for these two composers a more uniform distribution among the different correlation profiles, which also indicates a higher variability of composing schemes. For the other composers α 2 < α 1 is more recurrent. The appearance of two scaling exponents can be attributed to the presence of motifs as well as longer sections in the structure of the musical piece. Stronger short-range correlations (α 2 < α 1 ) are evidence of similar motifs repeated throughout the piece, while weaker shortrange correlations (α 2 > α 1 ) arise when the piece incorporates more varied short time patterns. In both cases, the structure given by the sections of the musical piece is reflected at long timescales. In a unique scaling case (α 2 = α 1 ), the musical patterns are preserved at all timescales, from motifs to phrases and sections. Even though there is evidence of correlations in their fluctuations, we are unable to extend this analysis to the remaining profiles (4 and 5) due to their scarcity.
To verify that these results are not a particularity of the multivariate application of the DFA method, we also applied DFA individually to each voice, finding quantitatively similar results for the various log(F(s)) functions (electronic supplementary material, figure S2). As mentioned before, scaling information can also be obtained from the series power spectra. In the supplementary material (electronic supplementary material, figure S3), as an example, we show cases of three different profiles. The localization of the crossover is consistent with the DFA calculations, however, in the latter they are better defined and more precise. The expected scaling relations among the exponents are corroborated.

Nonlinearity
Music, as well as language, are prototypes of so-called 'complex systems', which are frequently governed by nonlinear interrelations. Nonlinearity of a time series has been related to multifractality [26] and previous studies have already reported multifractal properties of music [5,28,29]. In view of these results, we search for nonlinear autocorrelations in the present work. To this end, we generate surrogate data and apply DFA to the 'magnitude series' [24] of both original and surrogate time series. The nonlinearity test is presented in the material and methods section.
To check whether the algorithm used for the generation of the surrogate data works well, we first compare the fluctuation functions derived from the original scores with those estimated from the corresponding surrogates. Owing to the fact that IAAFT-surrogates share the same linear univariate properties as the original times series, one expects that the results of the DFA from original pieces fall within the statistical range of the surrogate DFA results. In figure 11, we show examples of some of the results, with a 5% statistical significance level. One observes that all DFA estimations from the original musical pieces are within the area covered by the statistical fluctuations of the F(s) obtained for the surrogates. The same is true for all the remaining compositions treated in our study. This confirms that linear autocorrelations are conserved with high precision, since nonlinear properties will be destroyed by Fourier phase randomization. Hence, such surrogate data reflect appropriately the null hypothesis.
By inspecting the results obtained for the MDFA for some of the pieces in figure 12, one observes a marked difference between the surrogate data and the original series, viz. a clear evidence for the presence of nonlinear correlations in the time series of the music scores. Some of the fluctuation functions of the original data are outside the region covered by the surrogate data within a restricted range of timescales (like the compositions by Bach shown in figure 12); others show striking differences over the whole range of s (e.g. those of Beethoven).
The behaviour of the MDFA function of the original time series is slightly different from one composer to another. In Bach, the nonlinear correlations are more evident within shorter ranges; in Palestrina the behaviour is similar but the range of the nonlinear correlations seems somewhat longer than in Bach.

Discussion
In this study, we performed a statistical analysis of musical scores from several periods, focused on the behaviour of correlations, looking into their range, scaling properties and also a means for the detection of nonlinearity. To this end, we used DFA, which circumvents the question of whether musical pieces are stationary, which is important in short time series such as these. Furthermore, it provides a uniform approach for the detection and characterization of linear as well as nonlinear autocorrelations. With the implementation of a modified DFA for intervals and comparison with surrogate data, we tested the presence of nonlinear traits.

Regarding the power laws
We found that though a considerable fraction of all compositions treated in this study (38%) show a clear power law scaling over the whole range of possible timescales, the dominant profile (53%) shows a crossover between two scaling regimes. This observation is corroborated by the determination of DFA exponents. The fact that the strength of short-range autocorrelations might be different from that of longrange dependencies provides an indicator for qualitatively different composition structures. In our study, it turned out that usually short-range correlations are more pronounced than autocorrelations in long timescales. Although, via the estimation of Hurst exponents, we have so far been unable to determine a certain musical period or the style of a particular composer, some trends can be identified and for some composers specific traits have been uncovered. For example, in Shostakovich we find that pieces with stronger long-range correlations are almost as frequent as those with stronger short-range correlations; this feature seems to reflect a particular style of this composer. The above is even more pronounced for Dvorak, although with poorer statistics; here we found that pieces with stronger long-range correlations are twice as frequent as those with dominant short-range autocorrelations. Another clear example comes from Palestrina, for whom the finding of a marked dominance of short-range correlations over long ones appears to be a hallmark of his style. In general, the location of the distribution of the musical pieces within the (α 1 ; α 2 )-plane is revealing. Another interesting aspect in the context of a possible classifier is the presence of nonlinear autocorrelations. The fluctuation analysis of the magnitude series, together with an adequately designed surrogate test, provides an efficient tool for the detection and characterization of nonlinear dependencies within musical scores. Here, we find strong hints of peculiar, composer-specific, nonlinear features. surrogates. In the two functions of the piece from Bach, there is no significant difference between the original and surrogate data, while the difference in the MDFA from Beethoven is evident.
Significant differences in magnitude and temporal range with the results obtained from surrogate data may serve as indicators for specific temporal structures of a musical piece. Though these findings are not conclusive due the poor statistics, we strongly believe that linear as well as nonlinear autocorrelation structures are prominent candidates for characterizing certain epochs or, equivalently, particular composer styles. Further analysis in this direction is surely required.

The pleasantness of the nonlinear correlations
Since the work of Voss & Clarke [1] on music audio voltage, there have been numerous studies that refer to 1/f correlated noise in music as the most pleasant to the human ear [4,9,10,30]. In this case, the relation between fluctuations at different scales is preserved with a power law behaviour. There is a scale invariance in the power spectrum which by the Wiener-Khinchin relation is also manifest in the linear autocorrelation. On the other hand it has been argued by Schoenberg [31] that in order to optimize aesthetic appreciation of a musical piece, a certain equilibrium between regular time structures and variation should be encountered. Though scaling in linear autocorrelations may play a crucial role in this regard, regularity is not necessarily produced by self-similarity of different fragments of a musical piece. Nonlinear auto-interrelationships might generate less obvious, but possibly not less fascinating regular features. However, to the best of our knowledge, there are so far no studies mentioning nonlinear correlations in music. We designed an experiment in order to probe the aesthetic quality of the presence of nonlinear autocorrelations in compositions. To this end, we selected two different pieces and their corresponding surrogates ( figure 13). The MDFA of one composition created by Bach evidences only a weak contribution of nonlinear correlations, which are furthermore exclusively present on short timescales. The other opus from Beethoven shows nonlinear features on all scales, whose magnitude is growing with s. The selected pieces are the prelude No. 6 of the well-tempered clavier from Bach and the finale of the 13th string quartet from Beethoven.
We constructed MIDI sonifications of the two original pieces and their respective surrogates and conducted a survey for a quantitative evaluation of the aesthetic quality of each. In total, 1281 persons were consulted on how pleasant they perceived each of them in a scale from 1 (highly unpleasant) to 10 (most pleasant); for more details on the experiment see the electronic supplementary material. The results are summarized in figure 14.  In both cases, the evaluation of the original pieces turned out to be quite positive, with an accumulation of the scores between 8 and 10. However, the evaluation of the surrogate pieces was much more surprising. Here, the surrogate composition by Bach received significantly less low qualifications (1 and 2) and simultaneously more higher scores (in particular a score of 8) than in Beethoven.
The lack of nonlinear correlations in the original piece of Bach means that most of the regular structure is incorporated in the power spectrum. Self-similar musical motives are preserved under the conservation of the power spectrum. Therefore, also the surrogates maintain somehow the equilibrium between order and disorder of the original piece and encounter a benevolent evaluation.
The situation is different in the case of Beethoven's 13th string quartet. According to our results, this composition contains an important amount of nonlinear correlations, which are destroyed in the surrogate music. Under these conditions, the power spectra alone no longer assures the fine coordination between ordered repetition of motives and structural variability. Nonlinear correlations play now a major role and the musical piece differs much more from the noisy character of the surrogates. This last statement is true not only for the obvious case of Beethoven's 13th string quartet, where nonlinear features are quite pronounced, but also for the prelude No. 6 of the well-tempered clavier from Bach, where nonlinear correlations are seemingly present in a subtle manner; they play nevertheless an important role for the aesthetic appreciation. Also, in this case we measured a striking difference between the evaluation of the original piece and the surrogate replica. Hence, we may conclude that, in general, aesthetic perception in musical compositions depends importantly on their nonlinear autocorrelation time structures.

About algorithmic composition
One of the main interests for understanding musical structures is the development of algorithmic composition. Markov models, neural networks and generative grammars have been proved to be successful in generating note sequences with musical meaning at short timescales [32][33][34]. However, the research in the development of new models with short-and long-range correlations is still of interest in statistical physics [35,36]. The characterization of linear as well as nonlinear autocorrelations and the timescales where they are present could lead to more realistic models. The inclusion of nonlinear scaling behaviour by means of stochastic models could be helpful in the development of new techniques in algorithmic composition.

Conclusion
We have presented an analysis of music scores using DFA. We were able to identify different profiles in the fluctuation functions, which could be used for further classification of musical pieces. We found that not all the pieces have simple scale invariance in their fluctuations function; indeed the presence of a crossover between two different scaling regimes is more frequent. This evidences that pieces have different statistical behaviours within different ranges, i.e. different correlations at different timescales. By looking into the correlation profiles, we uncovered traits in the composition styles of some of the musicians. We exemplified the relation of structural modules with correlations at different timescales in two pieces. We also uncovered clear tendencies in composers and temporal evolution of compositions related with the cumulative distributions of the α exponents (in the case of single scaling) and with the mean and dispersion of the α 1 and α 2 exponents (in the crossover cases).
We further applied the MDFA method to the music scores and found evidence of the presence of nonlinear correlations. Similar to the linear DFA approach, different profiles for F(s)/s were encountered. We were unable, so far, to establish a relationship between the scaling profiles obtained by DFA and those of MDFA. We constructed surrogate pieces preserving the linear correlations of the original compositions and undertook a survey in order to evaluate the pleasantness of the surrogates. We found that nonlinear correlations could play an important role in the aesthetic appreciation of the musical pieces.
One of the aims of this paper was to contribute to establish criteria for the classification of musical compositions; some progress was achieved in this respect with the analysis of the different profiles identified in the DFA function. Additionally, our study provides elements and tools for the analysis of specific pieces. We believe this approach to the analysis of music scores, which unravels characteristics of composers, musical forms and periods, has the potential of contributing to the understanding of further issues such as musical evolution, composition and perception.
Data accessibility. The data and the code used in this work, as well as the electronic supplementary material are deposited at Dryad: https://doi.org/10.5061/dryad.6737v [37].