2.1 Challenges

Each NUScon challenge began with an unpublished experiment that was collected with US. The Challenge Committee was tasked with identifying appropriate experimental data that are high quality and representative of molecular sizes and experiment types of interest to our community. The descriptions provided to the contestants for each molecule used in a challenge data set are given in Table 1. A configuration file was deposited with each experiment, which provided basic metadata (e.g., number of points along the indirect dimensions, spectral width, and region of interest to be extracted along the direct dimension) that were used by the processing workflow. An example configuration file is included in the Supplement. We also include ¹H,¹⁵N projections from the 3D challenge data in the Supplement; these were not provided to contestants during the open challenge.

Table 1Descriptions provided to the contestants for each of the molecules used in the challenge data.

Download Print Version | Download XLSX

2.2 Synthetic peaks

Evaluating the performance of a spectral reconstruction method is a challenging proposition, especially if the procedure must be automated. The central challenge is that NMR spectra are used for many purposes, and the properties that make a reconstruction most fit for purpose are not always the same. Also, in practice, parameters of interest are generally extracted from spectra rather than from time-domain data. So, while it is possible to make an unbiased, automated evaluation of agreement between observed time-domain data and the corresponding time-domain version of a reconstruction, this does not necessarily capture how effectively a reconstruction may serve a particular spectral analysis task.

Along these lines, while it may be possible to measure fidelity of frequency, amplitude, and signal envelope (decay) of individual signals, this does not necessarily capture the fitness of a reconstruction to suit a given use. In particular, while NUS reconstruction methods can have excellent lineshape fidelity (Stern et al., 2002; Bostock and Nietlispach, 2017; Roginkin et al., 2020), it is seldom important in practice that a lineshape be conserved by a reconstruction, hence the common use of apodization that alters the lineshape, even if it conserves the integral (Naylor and Tahic, 2007). This is particularly true in multidimensional biomolecular data, where the indirect dimensions have little or no decay, so that the observed lineshape is primarily due to signal processing details rather than to any property of the underlying time-domain data. Furthermore, many applications, such as typical backbone assignment tasks, do not use peak heights or integrals quantitatively. Therefore, as a general point, it is often more useful to sharpen peaks or to suppress reconstruction artifacts than it is to faithfully reproduce lineshape.

The NUScon Planning Committee developed strong consensus that the single most important property of spectral reconstruction is for a real peak to be detectable, a fundamental requirement for all subsequent analyses. This, in turn, means that the evaluation of a reconstruction method is conflated with the procedure used to detect a peak. This is especially problematic for NMR, since the final decision about whether a spectral feature is a peak is often performed visually. Evaluation of the fidelity of the frequency, amplitude, and decay of a signal is likewise conflated with the procedure used to extract these parameters. The performance of such procedures is influenced by the presence and dynamic range of other peaks in the spectrum and by residual phase distortion, systematic artifacts, thermal noise, and the true functional form of the underlying NMR signals. Regarding this latter point, evaluations are further complicated by the potential for a single multidimensional peak to actually be a composite of mixed-phase signals in the presence of one or more unresolved couplings (Mulder et al., 2011) or the mixtures of signals with different relaxation pathways (Pervushin et al., 1997).

Keeping these challenges in mind, we developed a flexible procedure to automatically define and insert synthetic signals into existing time-domain data. These synthetic signals have a similar appearance and arrangement to signals in the original data. For example, signals in an HNCA spectrum will generally occur in pairs at the same ¹H,¹⁵N chemical shift and will share the same lineshape in the ¹H and ¹⁵N dimensions. Accordingly, the signal injection procedure provides options to specify how the synthetic signals should be placed relative to existing peaks. For example, we can specify that two synthetic HNCA peaks are inserted at a random ¹H,¹⁵N position where no peaks exist anywhere in the original 3D data but also that these two peaks are placed at ¹³C positions that correspond to peaks elsewhere in the original data. In order to better mimic the appearance of measured data and to help minimize any evaluation bias due to use of synthetic data, the individual signals allow for random perturbations of linewidth and phase and optionally one or more random unresolved couplings in each dimension. It should be noted that there are opportunities to improve the realism of simulations still further, including using composite signals of mixed phases or linewidths or by attempting to match the coupling constants and decays of synthetic peaks to the values of related peaks in the measured spectrum. For example, in the current HNCA peak generation procedure, peaks can be inserted to align with a ¹³C chemical shift of a peak elsewhere in the measured spectrum, and that synthetic peak will have randomly chosen ¹³C decay and couplings. In a more realistic simulation, the synthetic peak could be generated to match not only the ¹³C chemical shift of a peak elsewhere in the spectrum, but also that measured peak's ¹³C decay and couplings, provided that these details are known.

The signal injection procedure automatically characterizes the experimental reference spectrum to prepare for the injection of synthetic peaks. The procedure was built from C-shell and TCL scripts using the facilities of NMRPipe and its custom TCL interpreter nmrWish (Delaglio et al., 1995; Johnson and Blevins, 1994; Ousterhout, 1994). Here, we present an outline of an HNCA example used in the first NUScon, with complete data and scripts available on NMRbox (Maciejewski et al., 2017). The input to the procedure was the measured reference 3D time-domain data, which were used to establish acquisition parameters used for simulating time-domain data, and the corresponding 3D spectrum, which was used to establish the phase values for simulated signals so that the simulated data can be processed by the exact same scheme applied to the measured data. The output of the procedure was a table of synthetic 3D signals, a 3D time-domain simulation with only the synthetic signals, and a version of the reference 3D time-domain data with the synthetic signals added.

1. The fully sampled 3D time-domain reference data are Fourier processed to generate a corresponding spectrum.

2. Automated peak detection (detection of local maxima) was applied, and the largest peaks were retained (256 peaks for this HNCA example). The NMRPipe peak detection procedure estimated heights and linewidths by a multidimensional parabolic fit of the points around the maxima.

3. Using knowledge of median linewidths from peak detection, each point in the 3D spectrum was evaluated to decide whether it was in an “empty” region, meaning that there are no substantial signals within ±2 median linewidths. The criteria are that the rms of the region must be below 1.15 times the estimated noise and that no points in the region are above 5.0 times the estimated noise. The results were saved as a mask of the 3D spectrum, where intensities were set to 1 for each point with an empty neighborhood and 0 otherwise. Two-dimensional projections of the 3D mask were also prepared as shown in Fig. 2.

4. Using the peak and mask information from the steps above, the utility genSimTab.tcl was used to generate a synthetic peak table. Synthetic signals were only inserted in empty regions of the original data. Each row in the output peak table described one signal by a single amplitude and by a frequency, exponential decay, and phase for each dimension. Optionally, one or more couplings (cosine modulations) for each dimension may be specified. An example NMRPipe command illustrating the use of genSimTab.tcl is shown in Fig. 3. A key feature of genSimTab.tcl is that the desired properties of the signals to insert can be specified by a “program” of keywords and values, allowing substantial flexibility. Each program generates one or more “strips” (3D peaks at the same ¹H,¹⁵N coordinate). For HNCA data, each strip of synthetic signals consists of two 3D peaks. For example, the program

   $$\begin{array}{}\text{(1)}& \mathtt{count}=\mathtt{1},\mathtt{xy}=\mathtt{novel},\mathtt{z}=\mathtt{existing}\end{array}$$

   means generate one strip (count=1) of two peaks (HNCA peaks are generated in pairs) at a randomly chosen ¹H,¹⁵N location where no peaks exist anywhere in the original 3D data (xy=novel) but at ¹³C positions of known peaks elsewhere in the spectrum (z=existing).

5. Synthetic time-domain data were generated by the simTimeND application using the synthetic peak table and the acquisition parameters from the measured data as inputs. This included the placement of a group delay of “negative time” points at the start of each 1D time-domain vector if the reference data did not properly account for digital oversampling, as is the case for data from Bruker spectrometers (Moskau, 2002). The synthetic signals had exponential decays, with options for mixed-time data (Ying et al., 2007). The synthetic data also included a synthetic solvent signal of fixed width whose amplitude, phase, and frequency varied randomly in each 1D time-domain vector according to specified lower and upper bounds.

6. The synthetic time-domain data were added to the experimentally measured reference time-domain data. Example 2D projections and strip plots are shown in Figs. 2 and 4, respectively.

7. Given an appropriate sampling schedule, the fully sampled time-domain data can be resampled to generate NUS data via the nusCompress.tcl utility.

An advantage of this scheme is that it produces a spectrum processed identically to the experimental data but containing only synthetic peaks. This allows the simulated signals to be subjected to any desired detection or quantification method in either the presence or absence of signals and noise from the measured data.

Figure 2Synthetic peaks injected into HNCA empirical data. Two-dimensional projection of the 3D empty region mask (gray) overlayed with the 2D projections of the 3D reference spectrum (blue) and simulated signals (green). The areas shaded gray are regions in the measured data determined to be signal-free, such that a signal inserted at any shaded position will have no substantial signals nearby anywhere in the 3D spectrum. Synthetic signals were generated by the command in Fig. 3, which generates 10 strips of synthetic peaks at the ¹H,¹⁵N positions labeled here. Strips from positions 7 to 10 are shown in Fig. 4. As can be seen here, strips 1 to 7 are at 2D coordinates with no existing peaks, strip 9 is at a 2D coordinate that partially overlaps with existing peaks, and strips 8 and 10 are at 2D coordinates that correspond exactly to existing peaks.

Download

Figure 3Example genSimTab.tcl command. This command generates a table of synthetic peaks for an HNCA spectrum.

Download

Figure 4Strip plots of synthetic peaks. These strips show the auto-generated synthetic peaks and the corresponding measured data from the command in Fig. 3. The command uses five “programs” to produce strips of synthetic 3D signals at 10 ¹H,¹⁵N positions. The last four groups of strips (groups 7 to 10) are shown here, from the measured spectrum, simulated spectrum, and sum, labeled FT, SIM, and SUM, respectively. For these HNCA data, two simulated signals are inserted in every case. Strip group 7 corresponds to a program which inserts signals at ¹H,¹⁵N locations where there are no peaks anywhere in the original data, so that the “FT” strip in group 7 shows no signals. Strip group 8 corresponds to a program that inserts non-overlapping peaks at a ¹H,¹⁵N position where there are signals at other ¹³C positions in the original data. Strip group 9 corresponds to a program that inserts peaks which are offset by 0.5 ppm in the ¹⁵N dimension relative to an existing peak in the original data. Strip group 10 corresponds to a program which inserts two peaks that are 0.9 ppm apart in the ¹³C dimension. It can also be seen that the measured data for this case have a combination of peaks with a broad appearance, as in FT strip 8, and peaks with a narrow appearance, as in FT strips 9 and 10.

Download

2.3 NUS sample schedules

The NUScon software generates sample schedules by passing parameters to the nus-tool software package on NMRbox. The NUScon contest employed an exponentially biased sample schedule at a low and high coverage for each experiment. The parameters that defined each sample schedule are described in Table 2. The public archive of NUScon data contains the sample schedules along with the complete nus-tool command used to generate each schedule, including the random seed.

Table 2 has a column that shows the “grid” size, which is the dimensions of the Nyquist grid that spans the indirect dimensions. Following the notation established in Monajemi et al. (2017), each point in a sampling grid is an indel (indirect element). The column in Table 2 labeled “Indels” shows the number of indels that were collected by the corresponding sample schedule. Each indel represents four FIDs for these three-dimensional experiments with quadrature detection employed. All schedules used “full component sampling,” so all four FIDs were collected at each indel. In future work, we may consider partial component sampling (PCS; Schuyler et al., 2013), where the group of FIDs at each indel may be subsampled.

Table 2Challenge experiments and sample schedules. The left half of the table shows the NMR experiments that define each challenge data set. The right half of the table shows the sampling parameters used to generate a lower and a higher coverage exponentially biased sample schedule for each experiment. The column “Nuclei” lists the indirect dimensions, as specified in the fid.com scripts. The term “Indels” (i.e., “indirect elements”) is a single combination of evolution times taken across all indirect dimensions. Sample schedules customarily list the indels at which FIDs are collected. Dimensions that are constant time are indicated by zero-valued decay rates.

Download Print Version | Download XLSX

The nus-tool command to generate each sample schedule used input flags to force the collection of all FIDs at the “first” and “last” indels (i.e., the FIDs taken at the lowest and highest combinations of evolution times). The first indel, which contains the largest signal intensity, would very likely be sampled by any reasonable exponentially biased sampling schedule, but forcing its inclusion is an advantageous convention. The size of the Nyquist grid was explicitly defined in the NUScon workflow metadata. However, in practice, empirical data may not include such information, and the largest time increments observed in a sample schedule are often the best indicators of the size of the Nyquist grid. To comply with that convention, we forced the inclusion of the last indel.

2.4 Spectral reconstruction

We provided contestants with a template script for submission that accepts NUS time-domain data, a sample schedule, and various other parameters as inputs. Contestants were charged with using any of the tools on the NMRbox platform to process the data and produce a spectrum, subject to a basic set of rules. These rules include general guidelines, like the ethics statement “Any script that is not `reasonable' or tries to cheat the challenge problem is disqualified.” The rules also contain specific processing guidelines to ensure fair evaluations: “The spectrum produced by the user script must not exceed a final size larger than `rounding up' each dimension to the next Fourier number and 3 zero fills.” The rules document is available in the NUScon archive on NMRbox.

The range of reconstruction software programs employed is summarized in Table 3. Several of the software packages provide access to various reconstruction algorithms, and some contestants used different algorithms for different challenges. The specific methods and parameters used by each contestant for each challenge should be referenced directly in their submission scripts inside the NUScon archive. It is worth noting that the SEER (from Hansen) and low-rank (from Qu and Lu) software packages were new to NMRbox and were provided along with the contest submissions. This showcases the ability of the NUScon workflow to promote the exploration of novel techniques and the agility of the NMRbox platform to readily incorporate them.

Table 3Contestants and categorization of their chosen spectral reconstruction methods.

Download Print Version | Download XLSX

2.5 Peak picking

The peak picking step in the workflow was simple to deliver but deceptively complex to interpret. We used the peak picker from NMRPipe, but setting the threshold and controlling the number of peaks recovered required some attention. The recovered peak table was used by the metrics discussed in the next section and, as will become clear, the number of peaks in the recovered table had an effect on the resulting metrics. As a consequence, the threshold used when peak picking had a direct effect on the metrics. Given that NUS reconstruction methods are nonlinear, there is no single intensity threshold that can be used across all reconstructions for a given data set. Therefore, the peak picking was performed with a low threshold; the resulting peak table was sorted by peak intensity and truncated to include 1500 of the most intense peaks. The fixed number was determined by manual inspection of the uniformly sampled empirical data and was set conservatively high to ensure that all true peaks would be captured.

The current approach based on a fixed number of recovered peaks was effective but will be replaced by the use of the in situ receiver operator characteristic (IROC, Zambrello et al., 2017) method in future rounds. The IROC method is advantageous since it continuously varies the peak picking threshold from the most intense peak down to the noise and reports the recovery rate and false discovery rate at each threshold. The resulting set of data defines the quality of the reconstruction independently of a specific peak picking threshold.

2.6 Metrics

The NUScon metrics were designed to quantify how well reconstruction workflows recover synthetically injected peaks. The metrics are briefly described here, and their full mathematical definitions are presented in Appendix A.

• M1: frequency accuracy. A symmetric Hausdorff metric was used with a maximum distance penalty to determine the accuracy of the recovered chemical shift positions. This metric allowed for different numbers of peaks in the injected and recovered sets (i.e., it does not require a one-to-one correspondence between the sets).

• M2: linearity of peak intensity. This quantifies how well the intensities of the recovered peaks were mapped by a linear function to their corresponding injected intensities. The correlation coefficient was calculated from the NumPy python package.

• M3: true positive rate. This reports the percentage of injected peaks that were recovered.

• M4: false positive rate. This reports the percentage of recovered peaks that were false.

• M5: linearity of peak intensity. This metric delivers the same evaluation as M2 but was implemented using an NMRPipe function.

2.7 Rank lists

A rank list was generated for each metric applied to each combination of experiment data, synthetic peak table, and NUS schedule. These rank lists are in the NUScon archive on NMRbox in comma-separated value (CSV) files. The report function in the NUScon software package can be used to extract arbitrary subsets of these rank lists to assemble aggregate results, which is how the NUScon challenges were evaluated.