Piano multipitch estimation using sparse coding embedded deep learning

As the foundation of many applications, multipitch estimation problem has always been the focus of acoustic music processing; however, existing algorithms perform deficiently due to its complexity. In this paper, we employ deep learning to address piano multipitch estimation problem by proposing MPENet based on a novel multimodal sparse incoherent non-negative matrix factorization (NMF) layer. This layer originates from a multimodal NMF problem with Lorentzian-BlockFrobenius sparsity constraint and incoherentness regularization. Experiments show that MPENet achieves state-of-the-art performance (83.65% F-measure for polyphony level 6) on RAND subset of MAPS dataset. MPENet enables NMF to do online learning and accomplishes multi-label classification by using only monophonic samples as training data. In addition, our layer algorithms can be easily modified and redeveloped for a wide variety of problems.

3f 2f 7 ≈ 2, refer to the 9th row of Table 1). Analogically, we have the interval from f 3 's fourth octave to f 's 19th overtone is about − 2 cents (1200 * log 2 19f 2 4 f 3 ≈ − 2, refer to the 4th row of Table 1). One can easily testify the rest of the table.
Besides, the acoustical characteristics of different instruments make the problem even more difficult: on the one hand, timbre variation results in different overtone magnitude distributions; on the other hand, inharmonicity 2 leads to various overtone frequency distributions [11]. Pianos are especially harder to deal with than other stringed instruments due to the complicated way of strings being wired. Due to the inharmonicity and its uniqueness on different strings [11], the slight frequency mismatch between the first overtone of a note and its upper octave will cause an interference pattern (a.k.a. acoustic beat) if pianos are tuned by exact equal temperament. In order to eliminate such acoustic beats, pianos are usually tuned individually by well-trained experts (called Table 1 Frequency relationships between the overtones of a reference note n and its upper octave Numbers under "Harmonic" indicate the overtone indices of n. Column indices of "Harmonic" indicate the octave numbers starting from 0. Variance cents in the last column are rounded up into integers harmonic tuning, the deviation from the exact equal temperament often forms a Railsback curve [12]). The other challenge of MPE comes from the complexity of note combination. Strategies for solving multi-label classification can be generally categorized into two, "one vs. all" and "one vs. one, " respectively. Let the class number be n, the former needs n classifiers while the latter needs n 2 = n 2 −n 2 ones. Although it is computationally feasible for most circumstances, classifiers are trained independently from feature extraction. The lack of supervision in feature extraction may degrade the performance since it is more meaningful for features to minimize classification error rather than reconstruction error [13]. Another existing strategy needs 2 n classifiers by encoding multilabels into single-labels. It is only feasible when n is not large; otherwise, one may suffer from dimension explosion problem. Taking the piano for example, choosing 7 notes from 88 yields 88 7 ≈ 6.3 × 10 9 combinations. Even if only timbre and decay are included, it is almost impossible to construct and train such a large-scale dataset. Moreover, as one of the most commonly used features in acoustic music processing applications, time-frequency representation is constrained by the uncertainty principle. The algorithm performance then may be degraded by such deficient feature. Meanwhile, recent results have shown that feature fusion from different sensors (namely modality, one may consider someone's fingerprint and iris, or footages of some action from different angles) has advantages for recognition tasks (cf. [14][15][16] and references therein). Combined information from multiple sources is more robust and tolerant to noises and errors. Multimodal joint representation under constraints maximizes the utility of different features, which can be used more effectively in task-driven scenarios. Note that multimodal features are different from stacking multiple features into one because the latter does not take the modality relationship into account, and increasing dimensionality brings huge computation and storage costs.
Based on and inspired by the above discussion, we in this paper propose MPENet, which is a deep learning (DL, [17][18][19][20][21]) network enhanced by a novel multimodal sparse incoherent NMF layer (MSI-NMF layer). MPENet and MSI-NMF layer are implemented by Caffe [22]. Structures of training and test phase are given in Figs. 1 and 2, where tensors (i.e., Caffe blobs) are denoted by arrows, Caffe built-in layers are denoted by rectangles (computation) and ellipses (loss), layer collections by hexagons, our implemented layers by rounded rectangles. "Mod, " "cod, " and "pred" are abbreviations for modality, coding, and prediction, respectively. Modality indices are appended by dashes. Repeated elements (represented by dashed lines) are omitted for simplicity purpose. Network details are explained in Sections 3 and 4. MPENet incorporates the supervision from data and task to train dictionaries and classifiers adaptively and jointly. Representative and discriminative features then can be used to make multilabel inference directly from superposed inputs. Our main contributions include: • Lorentzian-BlockFrobenius sparsity: A novel · L−BF,γ is imposed to a multimodal NMF model. Penalty is determined by the magnitude of class templates of all modalities so that class sparsity can be ensured.

Related work
Owing to the non-negativity and superposition properties of musical spectra, Non-negative Matrix Factorization (NMF, [23]) is applied widely in the latest acoustic music processing studies. Musical spectral data is decomposed into a dictionary and corresponding coefficients (also referred to as codings, activations, or activities in some references, we may use any of them according to the context). Note that NMF algorithms converge in unsupervised fashion, only rank-1 decomposition makes sense for computational stability and uniqueness purpose. Thus, most methods utilizing NMF employ a three-step procedure: (1) training note templates individually from samples of each note, (2) constructing a dictionary by concatenating all note templates, and (3) estimating multiple notes by computing the codings with the dictionary fixed. In early studies, each template has only one atom (columns of a dictionary are called atoms). Weninger et al. [24] develop this simple structure by dividing note samples into two parts: onset and decay. Then, two atoms are learned respectively from both parts, which yields a two-column note template. Such dictionary helps to capture the feature variation and distinguish note state over time. O'Hanlon and Plumbley [25] take a further step on dictionary flexibility. Note templates are constructed by using linear combinations of several pre-defined fixed narrow-band harmonic atoms. The input spectral data is then approximated under β-divergence group sparsity constraint. Other methods employing similar idea but different implementations are proposed in [2,4,[26][27][28][29]. Such procedure uses fixed dictionary to get note activations during test, so MPE results heavily depend on the learned note templates, i.e., training samples. One has to retrain each template once new samples are added into training set. For other work using NMF with row/group sparsity and incoherent dictionaries, refer to [30][31][32] and references therein. Note that there are also studies that use unsupervised NMF instead of training note templates via isolated note samples. Bertin et al. [33] propose a tempering scheme favoring NMF with Itakura-Saito divergence to global minima. O'Hanlon and Sandler [34] propose an iterative hard thresholding approach for l 0 sparse NMF problem with Hellinger distance. ERBT spectrograms of polyphonic music pieces are decomposed directly and a pitch salience matrix is calculate to detect active notes. A semi-supervised NMF method can be referred to [35]. Many non-NMF based algorithms have been proposed for MPE problem. Tolonen and Karjalainen [36] divide the signal into two channels according to a fixed frequency and compute autocorrelation of the low channel and the envelope of the high channel to form summary autocorrelation function (SACF) and enhanced SACF (ESACF). The SACF and ESACF representations are used to observe the periodicities of the signal and estimate notes. Klapuri [37] calculates the salience representation through a weighted summation of overtone amplitudes. Three estimators based on direct, iterative, and joint strategies are proposed to extract notes from the salience function. Emiya et al. [1] employ a probabilistic spectral smoothness principle to iteratively estimate polyphonic content from a set of note candidates. An assumption of maximum number of concurrent notes (n max = 6) is imposed to avoid extracting overmany notes. Adalbjörnsson et al. [3] use a fixed dictionary to reconstruct input signal under block sparsity constraint. Notes are then identified through coding magnitudes. The fixed dictionary used here, however, is constructed according to equal-tempered scale so that the algorithm is unsuitable for instruments with inharmonicity.
Deep learning has been used to address AMT problem in recent papers. Sigtia et al. [38] presents a realtime model which introduces recurrent neural networks (RNN) into a convolutional neural network (CNN, with only convolution, pooling, and fully connected layers). Kelz et al. [39] compare the performances of networks with different types of inputs (spectrograms with linearly/logarithmically spaced bins, logarithmically scaled magnitude, and constant-Q transform), layers (dropout and batch normalization), and depths. Hawthorne et al. [40] propose a deep model with bidirectional long short term memory (BiLSTM) networks and two objective functions (onsets and frames), achieving state-of-the-art performance on MAPS [1] under configuration 2 described in [38]. For more acoustic music processing work using deep learning, refer to [40] and references therein. Note that the deep learning methods listed here all use music pieces as training data, which means polyphonic information can be accessed, hence music language model and classifiers are learned simultaneously.

Notation
Throughout this paper, we denote vectors and matrices by bold lowercase and uppercase letters, for example, v ∈ R m and M ∈ R m×n . Parts of vectors and matrices are denoted by subscripts: v i is the i-th entry of v; M i , M i→ , M i,j , and M i,j,p,q represent the i-th column, i-th row, (i, j)-th entry, and p × q block starting from (i, j)-th entry, respectively.

Prototype
In comparison with other information fusion techniques, multimodal joint sparse representation provides an efficient tool and results in superior performance [41]. Redundancy is generally employed in dictionary learning algorithms [15,42,43] so that training data can be fit better and codings can be more discriminative and sparser. Besides, l p,1 norm ( [15]) is usually used to regularize codings for row sparsity, where It enforces dictionaries of different modalities using same atom to present same event, for example, l 2,1 encourages collaboration among all modalities, and l 1,1 imposes extra sparsity within rows. For MPE problem, dictionary incoherentness should be imposed to provide flexibility of modeling universal note representations in contrast to redundancy. As we discussed in Section 1, single-atom note templates cannot cover the diversity of music spectra whereas NMF cannot guarantee the stability and uniqueness for multi-atom ones. Because harmonic tuning aggravates overlapping partials, we can not distinguish that a spectral peak is a note overtone or a summation of several ones, i.e., it is not feasible to decompose frequency domain into orthogonal bins according to the center frequencies of harmonic series. In order to detect notes directly from factorization, a "good" dictionary should be trained under the supervision of data and task, possessing the following properties: (1) note templates are mutually discriminative and (2) for a certain note, all possible variants can be and only can be represented by its templates.
Moreover, we improve l p,1 norm for two reasons. The first one is coding structure does not satisfy row-wise sparsity since dictionary incoherentness is imposed. Note samples are approximated by linear combination of its template atoms. The second one is l 1 norm imposes too much penalty so that every activation is either scaled down or zeroed out by soft threshold shrinkage [44]. For unknown number and loudness in MPE problem, each coding entry is crucial for detecting notes correctly, so we want to preserve as many effective activations as possible. Opposed to l 1 and l 2 , the Lorentzian-l 2 norm [45] penalizes large activations with small weights but the other way around so that non-zero activations keep their contributions. Besides, l 1 norm is not differentiable at 0, which makes the computation of gradient complicated ( [15] tackles this by introducing "active set"). The everywhere smoothness of Lorentzian-l 2 provides good convergence property. Figure 3 shows the contours of several common regularizations. Summing up the above discussion, the prototype of multimodal sparse incoherent NMF layer is a multimodal sparse incoherent NMF model whose cost function is, where superscripts indicate modality indices, m denotes modality number, f denotes feature dimensionality, n denotes class number, a denotes atom number of each class template, d = n × a is dictionary column number, and {μ, λ 1 , In (3), A i contains all template coefficients of the i-th class of all modalities. Frobenius norm incorporates the contributions of different modalities. Thus, Lorentzian-BlockFrobenius norm imposes class sparsity instead of row sparsity. The inner product term enforces that the dictionary columns have the least coherentness. It ensures the discrimination among class templates as well as the linear representation within each template. Note that analytic or straight optimization cannot be done for (3) because it is not jointly convex with respect to (w.r.t) D i , i = N m and A. It is convex w.r.t either one while the other fixed. Hence, many alternating schemes ( [42,44,46,47]) split (3) into two subproblems, sparse coding and dictionary learning, respectively.

Structure
MSI-NMF layer is constructed by re-translating the two subproblems of (3), where we treat A as layer outputs and D i , i = N m as layer parameters. Using the same network notations as in Section 1, Fig. 4 shows the structure of MSI-NMF layer, where layer parameters are denoted by texts within parentheses.

Forward pass
The forward pass produces the solution of a multimodal non-negative sparse coding problem whose cost function is defined as Then, (4) can be rewritten as where (5) can be solved using Alternating Direction Method of Multipliers (ADMM, ref. [44,46,47]), details are given in Algorithm 1 (proofs in Appendix), where is given in Algorithm 2.
Algorithm 1 Forward Pass of Multimodal Sparse Incoherent NMF Layer: Multimodal Non-negative Sparse Coding : until a converges Ensure: a

Backward pass
The backward pass is to update D through gradient descent. The incoherentness constraint is treated as weight decay of layer parameters. Denoting the network cost function by l net , the new cost function becomes Then, we have where p = N f , q = N md . According to the chain rule, the first term of (9) is ∂l net In order to get ∂a ∂D p,q , recalling that a is a minimizer of (5), taking the derivative w.r.t a, we have whereW = (a) is defined as

Algorithm 2
: Inner Update Algorithm of a in Algorithm 1 if λ 4γ 2 then 5: if γ 2 = 1 27 and λ = 4γ 2 then end if 8: else 9: if > 0 then goto 6 10: end if 13: end if 14: end if 15: end for 16: for j = 1, 2, . . . , md do a j = u j p k , k = j/a mod n 17: end for Ensure: aW a is defined in (7). Then, ∂l new ∂D p,q can be computed because ∂a ∂D p,q can be obtained by taking the derivative w.r.t D p,q on (11) and ∂l net ∂a is given by the last layer. The backward algorithm of proposed layer is listed in 5: 6: 7: (20) 8: 9: end for Ensure:

MPENet
In this section, we detailedly explain the layer and tensor specifics of MPENet. In order to avoid misunderstandings caused by layer names in different deep learning frameworks (for example, commonly called "fully connected" is named as "inner product" in Caffe, "linear" in PyTorch, and "dense" in TensorFlow), during illustration, we will give the mathematics expression of some layers if necessary. Meanwhile, in order to give the most direct ideas of how MPENet is constructed, we switch to Caffe terminology accordingly (see Figs. 1 and 2).
In Figs. 1 and 2, training and test phases have same core modules, differences only locate in the top layers. "Data" layer produces multimodal features and their labels. Labels are binary vectors whose entries are 1 if corresponding classes are active and 0 otherwise. "Sig-moidCrossEntropyLoss" layer is a stack of "sigmoid" layer and "cross-entropy" layer. Cross-entropy loss is defined as where p andp are predictions and labels.

Deep multi-label prediction module (DMP)
The structure of "Deep Multi-label Prediction" (DMP) is shown in Fig. 5. "Slicing" layer segments an n × a vector into n parts. "Detection" is a classifier module with replaceable structure, and is supposed to output the existence magnitude according to the input. "Concat" (concatenation) layer joints n detections to form a multilabel prediction. In our experiment, five layers are used to implement "Detection" module (structure is shown in Fig. 6). "InnerProduct" represents the transform from x ∈ R m to y ∈ R n y = Wx + b, W ∈ R n×m , b ∈ R n "ReLU" (Rectified Linear Unit) stands for the transform from x ∈ R m to y ∈ R m For other tasks, one can modify this combination accordingly.
The reason for such structure roots from the property of incoherent dictionaries. If samples of certain class can be and only can be represented by its template atoms, the existence of this class is only related to the coefficient magnitudes. "One vs. all" strategy can be employed natively. If dictionaries are not as good as expected, crossentropy loss will correct each "detection" module as well as dictionaries of the proposed layer through backward pass, which completes a positive circle.

Multi-label accuracy layer
Multi-label accuracy layer is implemented to conduct training and test in a unified framework. It consists of three sequential operations: sigmoid activation, binary output, and metric computation. The second one outputs either 0 or 1 according to the comparison result between the sigmoid activation and a predefined threshold t. The third one calculates the Precision (P), Recall (R), and F-measure (F) according to the binary outputs and the ground truths, where TP, FP, and FN stand for true positive, false positive, and false negative, respectively.

Experiment results
In this section, we first briefly demonstrate the dataset and features used in our experiment, then illustrate parameter initialization and network configuration in detail. Piano MPE results, experiment results about how MPENet works, timbre robustness results, and AMT results are given in the end of this section.

Dataset and features
MAPS [1] is a commonly used piano dataset for multipitch estimation and automatic transcription. It contains nine kinds of recording conditions (referred to as "StbgTGd2, " "AkPnBsdf, " "AkPnBcht, " "AkPnCGdD, " "AkPnStgb, " "Sptk-BGAm, " "SptkBGCl, " "ENSTDkAm, " and "ENSTDkCl"), two of them ("ENSTDkAm" and "ENSTDkCl") are from real pianos and seven are synthesized by softwares. Each kind has same subset hierarchies which include ISOL (monophonic recordings), RAND (random combination), and UCHO (chords). ISOL/NO subset, which contains 264 monophonic wav files covering 88 notes (n = 88) and 3 loudness levels, is used as training set. RAND subset, which contains 6 polyphony levels ranging from 2 to 7 (labeled as P2-P7), is used as test set. Each one of P2-P7 has 50 files, and the note combination of each file is generated randomly. In [1], a 93-ms frame which is 10 ms after onset of each file in P2-P6 is analyzed. As comparison, we conduct similar evaluation in our experiment. P7 is used as validation set for parameter tuning.
Each wav file in MAPS is stereo with sampling rate 44100 Hz. To extract features, we firstly generate a mono counterpart by averaging both channels. Then, the silent part of each counterpart is truncated according to the provided onsets and offsets. Finally, two kinds of features (m = 2) are extracted from the remainder by Short Time Fourier Transform (STFT) and Constant-Q Transform (CQT). The reason for using STFT and CQT is mainly because the former has good resolution in high-frequency domain while the latter does well in lowfrequency domain. STFT and CQT features are further transformed into non-negative dB scale using h(·) log 10 (·/ + 1) log 10 (1/ + 1) where · is either CQT or STFT feature and is the machine precision. Other extraction specifics are listed in Table 2, where flen, slen, minf, maxf, dim, ppo, and nfft are abbreviations for frame length, step length, minimal frequency, maximal frequency, dimensionality, partition per octave, and n-point Fast Fourier Transform, respectively.

Parameter initialization
Due to the non-convexity of problem (3), only local minimization can be guaranteed. The initial value of dictionaries is crucial for convergence and performance. Since totally random initialization makes the codings of first several epochs meaningless, it is a waste of time and computation resources. Plus, because each monophonic file in the training set lasts for over 2 s, many samples are similar to each other during decay. It is not reasonable to initialize the dictionary using random samples as most dictionary learning algorithms do [15] either. In our experiment, two procedures are employed to initialize the dictionary.
To avoid the heavy overhead caused by joint learning of dictionaries and classifiers, before really getting into MPENet, we propose a pre-learning phase called Label Consistent Incoherent Dictionary Learning (LCIDL) derived from [48,49] to obtain a better start than random initialization and sample initialization for dictionaries in MSI-NMF layer. The cost function of LCIDL is It is worth emphasizing that (23) is a plain data-driven problem, and neither MPENet nor classifiers are involved  at the time. It has nothing to do with deep learning and can be implemented by any language. The form of L is the extension of binary labels to impose classification information, because there are no note probabilities but only codings on our hands. Likewise, LCIDL also needs a good dictionary to start for acceleration. In order to find it and determine the atom number, a fast clustering algorithm based on density peaks [50] is employed to filtrate samples hierarchically. Specifically speaking, we first extract 30 cluster centers from each modality of each file in the training set. Then, we stack them according to their note indices. This gives us two matrices with 810 columns for each note (810 = 30 × 3 (loudness) ×9 (recording)). Finally, through computing density peaks on these two matrices and considering the overhead and efficiency of computation and storage, we empirically set the atom number of each note template to  (3)) and obtain a 576 × 1320 matrix and 648 × 1320 matrix for starting LCIDL. After LCIDL is done, we obtain a "roughly good" dictionary, it has low reconstruction loss, incoherentness, and coding shape like L as (23) governs. When the real training of MPENet begins, this "roughly good" dictionary is copied into MSI-NMF layer, and classifiers are initialized randomly. During the first several epochs of training, the learning rate of classifiers is relatively larger than that of MSI-NMF since we want to hold codings a little bit to fit classifiers first. As the classification error decreases, the learning rates of all layers become equal to do joint learning.

Network configuration
Choices of parameters used in MPENet are all empirical. The output numbers of three "InnerProduct" layers  in Fig. 6 are set to be 60, 30, and 1 from left to right. All three layers use bias term. For MSI-NMF layer, we use λ 1 = 0.15, λ 2 = 0.1, μ = 1.32, ρ = 0.2, γ = 1.09, and t = 0.5 for training. Considering that only one note is active at a time during training whereas at least two are active concurrently during test, test constraints should be weaker than those of training. Limited by the computation overhead, a fully greedy search cannot be done to get the best result. Therefore, we initialize several groups of parameters, and the one with λ 1 = 0.03, λ 2 = 0.1, μ = 1.32, ρ = 0.2, and γ = 0.55 gets the best result through evaluating the test net on P7. To tune binary threshold t of multi-label accuracy layer, we plot a Precision-Recall Curve in Fig. 7 according to the evaluation result on P7, where t ranges from 0 to 1 with step 0.001. Through the figure, we find that "detection" modules produce very polarized outputs. The best result (91.37% Precision and 70.92% Recall, see the gray square in Fig. 7) and similar ones can be achieved when t is within a relatively large interval around 0.5. Therefore, we keep t = 0.5 unchanged for test.

MPE results
Evaluation metrics are listed in Tables 3, 4, and 5 when the polyphony level is unknown. Our network outperforms all other algorithms on Recall and F-measure. Precisions of P2, P3, and P6 get the second best results with slight gaps compared to Emiya's. The decrease of sparsity constraint is the reason for this result shortage. During evaluation, we can achieve over 99.9% F-measure for P2 and P3 if we use training parameters. Such configuration can also maintain high Precision results for P4-P6; however, Recall will drop dramatically due to strong sparsity. It is a trade off and contradiction between sparsity and concurrent notes.
We also report the evaluation results in Table 6 when polyphony level is known as prior. For polyphony level k, we choose the indices of first k largest outputs in sigmoid layer as active notes. Results show that F-measure increases for P2 and P3 while things are different for P4-P6. It states a fact that when concurrent number is small, the ground truths have higher probabilities than others in our algorithm; as concurrent number grows, undetected ground truths become undetectable.

How MPENet works
In order to show how each key part of MPENet contributes to the performance, we conduct four groups and seven in total experiments. Considering the combination complexity, each experiment only changes single part to show its impact on the system. Group indices, names, and settings are listed in Table 7, where √ and × indicate presence and absence, respectively; the setting described in Sections 5.3 and 5.4 is called MPENet-default (MPENetd for short). Unless otherwise specified, all experiments in this subsection share same parameters with MPENetd except the modified part. Implementation details and results are explained in the following subsections. Table 7 as training inputs to show our multimodal efficacy. Results are plotted in Fig. 8. We find that STFT's Precision outperforms CQT's on all test sets, while the STFT's Recall decreases substantially from P3. MPENet-d, as expected, incorporates the advantages of both modalities and amend their drawbacks.

Atom number
Group 2 (Config.3-Config.5), in conjunction with MPENet-d, shows the influence of atom number on our system. We only change a described in Section 5.2 to initialize dictionaries with different sizes. Results are plotted in Fig. 9. Interestingly, the Precision of each one in group 2 gets improvement except P2 of Config.4. Especially, the    (5)) and the F-measure of P6, we consider MPENet-d sufficient enough. For those with unlimited computation resources, one can modify a and re-validate corresponding parameters for better performances.

Joint learning
Config.6 (without joint learning) is implemented by dividing dictionary and classifier learning as separate operations. During training, we first learn dictionaries by using Final metrics are the average over all pieces Results with asterisks are reimplemented by [40] (23) in Section 5.2 since it is the only way to impose supervision. Then, we compute the codings of monophonic training data by using (4) with learned dictionaries. Finally, monophonic codings are directly fed into "dmp" modules plotted in Fig. 1 to train classifiers. During test, we first compute the codings of polyphonic test data in P2-P6, then use the test phase in Fig. 2 to get metrics. Note that dictionaries are only learned once and do not change any more during coding computation and classifier training. Comparison results are shown in Fig. 10, where we find that although Config.6 beats MPENet-d by 2% constantly on Precision, the Recall of Config.6 has increasing gap compared with MPENet-d's as polyphony level grows. As a result, MPENet-d outperforms Config.6 greatly on F-measure for P4-P6.

Dictionary incoherentness
For Config.7 (without dictionary incoherentness), we remove the incoherentness regularization in (3). Due to the absence of incoherentness, block sparsity makes no sense then. The cost function used in Config.7 becomes where Lorentzian-Row_l 2 is defined as The corresponding form of (23) then becomes Forward and backward algorithm for (24) can be derived according to Algorithms 1 and 3. During training, we find the loss of training phase stays to a relatively high value (about two order higher than that of MPENet-d). Things do not change even if we reinitialize the parameters or train for extra several epochs. Moreover, during test, the sigmoid outputs of multi-label accuracy layer for P2-P6 are all less than the detection threshold t, so the metrics of Config.7 are all zero, test fails.
Summing up the results, we find that incoherentness regularization is crucial for MPENet while modality, atom number and joint learning only affect performance. If sorting them by importance, we have incoherentness modality joint learning atom number

Timbre robustness
In order to explore the generalization error of MPENet, another experiment is conducted by evaluating timbre robustness. In this experiment, we only choose "ENSTDkCl" as training set and test on other eight kinds.
Parameter settings are all kept the same as in the last subsection. The reason for choosing "ENSTDkCl" is that it is recorded from a real piano with "close" recording condition. Inharmonicity, tuning, timbre, decay, background noise and all other factors that can influence spectra may be very different from the other eight. The results in Figs. 11, 12, and 13 show that MPENet becomes overfitting since all metrics of test sets drop fairly except "ENSTDkAm" (only recording condition is different from training set).

AMT results
Due to the underlying strong relationship between MPE and AMT, and in order to further explore the capacity of MPENet, we also conduct AMT experiments following configuration 2 described in [38]. Specifically speaking, we use total 60 full-length music pieces contained in the "MUS" subsets of "ENSTDkCl" and "ENSTDkAm" as input and run MPENet frame by frame. Parameters and hyper-parameters are all kept the same as "MPENet-d" (c.f Sections 5.3 and 5.4). In line with the training phase of MPENet, the ground truths of music pieces are generated by discretizing note durations provided in corresponding txt files. Table 8 gives the frame-level average AMT results of MPENet, with comparison to state-ofthe-art performance reported in [40]. MPENet maintains the Precision as in Section 5.4, but performs poorly on Recall. Figures 14 and 15 reveal some occasions where false positives and false negatives take place, in which green indicates true positives, red indicates false negatives and blue indicates false positives. In brief, false positives consist of a few wrong chord detections and many scattered unrelated notes; while false negatives come from massive so-called super-combinations (number of simultaneously active notes is over 7) and plenty of notes with long duration. The former circumstance of false negatives, as discussed in Section 5.4, is an inevitable result caused by sparsity constraints. For the latter circumstance of false negatives, however, we think the reasons behind such behaviors are mainly caused by two aspects: (1) the training set of MPENet lacks negative samples. Since training loss is not zero and recall the decay properties of piano notes, classification errors of monophonic samples in training set include wrong detections and missing detections. The lack of negative samples makes MPENet tend to distinguish the beginning from the end of same note, which leads to insufficient durations; (2) MPENet knows nothing about music language, which prevents MPENet from rejecting scattered, unreasonable detections.
In order to compensate the second aspect discussed above and incorporate music prior with MPENet, we also train a simple recurrent network whose structure is given in Fig. 16. The training set follows the configuration 2 in [38], which consists of all music pieces in MAPS except the ones in "ENSTDkCl" and "ENSTDkAm". The music context is constructed as follows: for a certain time step t, the recurrent network takes current binary label l t as ground truth, and 16-time steps (about 100 ms) before and after it ({l t−16 , . . . , l t−1 , l t+1 , . . . , l t+16 }) as input. The number of hidden state in Gated Recurrent Units (GRU) is 200. The output number of all InnerProduct layers is 200 as well except the last one is 88 for label consistency. Before concatenation, we only keep the last step's output of each GRU because they accumulate all the information through time. After the training is done, we use this recurrent network to perform Gibbs sampling on the sigmoid output of MPENet. If letting the total frame number of test dataset be N and running Gibbs sampling N times as one full step, we find that one full step gives the best performance improvement (1.81% on Precision and 0.11% on Recall). As expected (see Fig. 17 for example), the recurrent network smooths out some scattered detections. With more than one full step, however, the recurrent network breaks down the initial MPENet output and tends to "recompose" it into a new piece of music, which results in a major metrics decreasing. Note that indices of Gibbs sampling are selected randomly, so not all frames have been updated during one full step. Also note that our recurrent network has way shorter memory than those in [38,40], so it learns little music language and only has effects on scatted detections. Because AMT is not the concern of this paper, we do not experiment more here but maybe focus on possible AMT-related refinements in the future work.

Conclusions
In this paper, we propose a new deep learning layer based on a NMF model with multimodal inputs under sparsity and incoherentness constraints. Such "layerization" of optimization problem provides the possibility to learn dictionaries and other features jointly under a unified deep learning framework. It enables modularization, online learning, and parameter fine-tuning for dictionary learning problem, which can be used to simplify the model refactoring and extension. In comparison with the "high level" features produced by other deep learning layers, the proposed layer learns discriminative and representative dictionaries so that the outputs are more realistically meaningful. Experiment results demonstrate that our test net improves the MPE performance substantially on MAPS dataset. Restricted by hardwares, we pay more attention to layer algorithm and the network structure than model training. Unlike those fully explored and well-tuned deep learning models, MPENet with empirical parameters, simple layer combinations and shallow structures have plenty room for improvements. For future work, there are several directions that can be considered: (1) from the layer point of view, performance grows with the increasing modality number. According to our experiment results, automatic parameter adaptation will also improve the estimation greatly; (2) from the network point of view, regularization, depth, and structure are new focuses for extracting more representative and robust features.

Endnotes
1 Scientific Pitch Notation is used to represent notes, i.e., sub-contra octave is indexed by 0.
2 For certain stringed instruments, overtones are close to but not exactly integer multiples of the fundamental frequency, the degree of departure from whole multiples is called inharmonicity.
⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎩ a = min a λ 1 a L−bl 2 ,γ + ρ 2 t − a + b 2 2 + δ R md 0 (a) (28) Solving (27) yields which is the update of t in Algorithm 1. To solve (28), we first change its form into a = min a λ 1 ρ a L−bl 2 ,γ + 1 2 t − a + b 2 2 + δ R md 0 (a) (30) Denoting λ = λ 1 ρ and using Karush-Kuhn-Tucker conditions, we introduce v ∈ R md 0 . The Lagrange function of (30) is and KKT conditions are where i = N md andW is defined in (13). It is easy to find (32) can be split into n independent groups and t i , b i , v i are defined accordingly. For any i = N n , we omit the subscript and let a = a i and u = t i + b i + v i , we have equations Substituting y 1 and y 2 into (47), we have where Note that ax 3 + bx 2 + cx + d = 0 having three different real roots, so its derivative 3a 2 + 2bx + c has two different real roots, i.e., 4b 2 − 12ac = 4A > 0 constantly. Finally, for θ ∈[ 0, π], it is easy to find that x 2 < x 3 < x 1 , we have β 1 = x 2 , β 3 = x 1 . Substituting the coefficients of (42) into x 2 and x 1 , one can get the equivalent expression described in Algorithm 2.
Back to v, according to (35), we have Combining the constraints of (34) and (33), we have v i = 0, Summing all the above discussion up completes Algorithm 1.

Proposition 2 ∂l new
∂D i , i = N m described in Algorithm 3 is the gradient of l new w.r.t D i .
Proof This proposition exploits the fact that the coding and dictionary of any two different modals are independent. First of all, (11) can be rewritten as equations Taking the derivative w.r.t D k i,j , we have where i = N f k , j = N d and E k ij ∈ R f k ×d denotes an all-zero matrix except the (i,j)-th element is 1.
Recalling the definition ofã in (7), then the q-th value of where ζ q = 2 γ 2 +ã q/a , q = N d Combing (52) and (53) and omitting some reduction and rearrangement, we have where V k is defined in (16). Let P k D kT D k +λ 1 W I−WV k V kT +λ 2 I, Q k P k −T ∂l net ∂A k According to (10), The first term of (55) is The second term is Substituting (56) and (57) into (55), we have where U k is defined in (17).