Frank Delaglio, Ph.D.

19804 Maycrest Way
Germantown MD 20876 USA

301 806-0867
delaglio@nmrscience.com



NMRPipe Processing Functions
MEMND: ND Maximum Entropy Calculation.

Flag Argument Default Description
 -ndim nDim 1 Number of Dimensions.
 -sigma sigma 1.0 Std Dev of Time-Domain Noise. Use Zero for Auto Estimates.
 -zf zfCount 2 Zero Fill Count. Time-Domain Input Only. (also use -xzf -yzf etc).
Reconstruction Mode Flags:
 -neg Use Two Channel (+/-) MEM (Default).
 -pos Use Positive-Only MEM.
 -zero Correct Zero Order Offset (Default).
 -nozero Suppress Offset Correction.
 -freq Output Freq-Domain Result (Default).
 -time Output Time-Domain Result.
Convolution Parameters (also use -x.. -y.. etc). Defines window used on MEM data:
 -conv {FILE fName} Keyword FILE followed by Convolution Function File Name.
 -conv cFn Convolution Window Function: SP EM GM GMB TM TRI JMOD.

 -cQ1 cQ1 Convolution APOD Function Parameter Q1.
 -cQ2 cQ2 Convolution APOD Function Parameter Q2.
 -cQ3 cQ3 Convolution APOD Function Parameter Q3.
 -csin Sets Q1 to 0.0 for Sine JMOD.
 -ccos Sets Q1 to 0.5 for Cosine JMOD.
Deconvolution Parameters (also use -x.. -y.. etc). Defines inverse window used on original data:
 -deco {FILE fName} Keyword FILE followed by Deconvolution Function File Name.
 -deco dFn Deconvolution Window Function: SP EM GM GMB TM TRI JMOD.
 -deco dFn None Deconvolution Window Function:
 -dQ1 dQ1 Deconvolution APOD Function Parameter Q1.
 -dQ2 dQ2 Deconvolution APOD Function Parameter Q2.
 -dQ3 dQ3 Deconvolution APOD Function Parameter Q3.
 -dsin Sets Q1 to 0.0 for Sine JMOD.
 -dcos Sets Q1 to 0.5 for Cosine JMOD.
Data Clipping and Scaling:
 -ts tsMode 0 Time-Domain Scaling Mode: 0=None 1=Divergence 2=Always.
 -clip cMode 0 Freq-Domain Clipping Mode: 0=None 1=Pos 2=Neg 3=Both.
 -min cMin -1e17 Minimum for Clipping.
 -max cMax +1e17 Maximum for Clipping.
MEM Parameters:
 -alpha alpha 1e-4 Scale for MEM Update.
 -tScale tScale 0.5 Scale for Amplification Threshold.
 -tLimit tLimit 1.5 Chi Limit for Use of Threshold.
 -cScale cScale 1e2 Desired Sigma After Scaling.
 -eScale eScale 1.0 Scale (A) for Entropy Calc.
 -lamb lamda 1e-5 Initial/Minimum Lambda Value.
 -step step 0.02 Scale for Lamda Update.
 -x0 x0 5.0 Second-Order MEM Threshold.
Convergence Tests:
 -cLim cLimit 1.0 Chi-Square/N Limit for Convergence. (same as old option -zFact)
 -zLim zLimit 1.0 Z-Statistic Limit For Convergence.
 -itmax itMax 50 Maximum Iteration Count.
 -itmin itMin 3 Minimum Iteration Count.
 -div divCount 3 Maximum Divergence Count.
 -report reportLev 0 Report Level: 0=No Report. 1=Final Iteration. 2=Each Iteration. 3=Save F and dC.
Notes:
  1. Important options are -ndim -sigma and -alpha.
  2. Underestimate time noise value for -sigma to avoid missing peaks.
  3. Use a minimum of 2 zero fills (size*4).
  4. Use -report 2 for statistics of each iterate.
  5. Reduce -alpha value if results diverge.
  6. Chi/N statistic should decrease to 1.0 in a good reconstruction.
  7. MEM Z statistic should increase to 1.0 in a good reconstruction, but values as small as 0.1 are not unusual.

MEM performs Maximum entropy reconstructions of 1D, 2D, and 3D spectra and spectral subsets (such as 2D planes of a 3D spectrum). The MEM reconstruction can be thought of as an attempt to create a spectrum of optimal quality or sharpness which is nevertheless consistent with the original time-domain data. The MEM function in nmrPipe is a two-channel implementation of Gull and Daniell's algorithm, with modifications introduced by Wu. The method and its basis is described in the following references:

Sibisi, S. (1983) Nature, 301, 134-136.
Gull, S.F. and Daniell, G.J. (1978) Nature, 272, 686-690.
Wu, N.L. (1984) Astron. Astrophys., 139, 555-557.
Laue, E.D., Skilling, J., Staunton, J., Sibisi, S. and Brereton, R. (1985) J. Magn Reson., 62, 437-452.
Hoch, J.C., Stern, A.S., Donoho, D.L. and Johnstone, I.M (1990) J. Magn. Reson., 86, 236-246.
Stephenson, M. (1988) Prog. NMR. Spectrosc., 20, 515-626.

USING MEM A typical MEM application might involve using 2D MEM on the planes from the indirect dimensions of a 3D spectrum. In this case, steps like the following will be performed:

  1. Decide if MEM is suitable for the application. It may be useful for improving resolution or reducing truncation artifacts; this particular implementation may not be suitable for data where precise quantification of peak heights or lineshapes is required, or in cases of high dynamic range. MEM may also be unstable when used on large data arrays, or simply impractical due to memory and computational requirements.
  2. Process the data completely with an ordinary Fourier transform scheme first, to determine details of phase correction, baseline adjustment, etc. Construct the scheme so that the dimensions of the result can be inverse-processed, and so that the axis-order of the result allows access to the data planes where MEM will be applied.
  3. Use nmrDraw to extract a plane from the indirect dimensions, and to estimate the noise level of this plane in the frequency-domain. If deconvolution will be used, estimate the intrinsic linewidths (linewidths of the data without window functions) as well.
  4. Use the program showApod to find the time-domain noise level based on the value of the frequency-domain noise level.
  5. Create and test a trial MEM processing scheme for the extracted 2D plane. Try executing the scheme without the MEM function first, to make sure that the details of inverse processing and reprocessing are correct.
  6. Execute and test the trial MEM scheme, adjusting MEM parameters as needed so that convergence is achieved in a reasonable number of iterations, and the MEM result seems to include all of the signals observed in the original data.
  7. Create and execute a similar MEM processing scheme for all the indirect dimension 2D planes in the 3D data, using the steps and parameters from the trial 2D scheme.

THE MEM ALGORITHM

The MEM algorithm implemented in nmrPipe is iterative. It starts with a trial MEM spectrum, which is initially uniform. At each iteration, the trial MEM spectrum is inverse Fourier-transformed and then subtracted from the original time-domain data to form a time-domain residual. The time-domain residual is then forward Fourier-transformed, and used to update the MEM spectrum. In the case of MEM deconvolution, the same procedure is used, except that the MEM data is broadened by a window function at each iteration, in order to generate a sharp reconstruction which when broadened is consistent with the original time-domain data.

A rough summary of the Gull and Daniell algorithm (without accommodations for combinations of positive and negative signals) follows.

D: is the original time-domain data.

n: is the total number of points in D.

sigma: is the standard deviation of the noise in D.

F: is the current MEM spectrum at a given iteration.

R: is the inverse Fourier transform of MEM spectrum F.

C: is the time-domain residual between the original data and the MEM data:

   C = (D - R)/(sigma^2)
U: is the forward Fourier transform of C. Therefore U can be thought of as the frequency-domain residual between the original data and the MEM data.

s: is the entropy value of MEM spectrum F, given by:

   s = -Sum( F ln(F/g) )
where g is a suitable scaling factor; s can be thought of as a measure of spectral quality in the sense of sharpness or simplicity.

chi2: is the chi-square value of the residual between the original data D and MEM data R:

   chi2 = Sum( (D - R)^2 / (sigma^2) )

The chi2 value should approach n when the MEM data is consistent with the original spectrum to within the noise level.

The MEM procedure attempts to find a spectrum F which will maximize the expression:

        q = s - lambda*chi2

where lambda is a La Grange multiplier which is chosen so that chi2 approaches n at the solution, or in other words so that the maximum entropy spectrum F is consistent with the original time-domain data D to within the noise. This leads to a solution spectrum B with scaling factor a where:

        B = a*exp( -1 + 2*lambda*U )

which corresponds to an exponential amplification of the frequency-domain residual U. In practice, this exponential form is unstable, so the solution spectrum B is found iteratively:

  1. Start with D (the original time-domain data with noise sigma), and F (the current MEM spectrum), which is initially uniform.
  2. Find R, the inverse Fourier transform of F, after using a Hilbert transform to reconstruct the imaginary part of F. If desired, apply a line broadening or other window to R in order to mimic the signal envelope in the original time-domain data D.
  3. Find the time-domain residual C, (D - R)/(sigma^2), as well as the chi-square value of the residual, chi2.
  4. Find U, the frequency-domain residual, by forward Fourier-transform of residual C.
  5. Find a suitable lambda value for the current iterate.
  6. Find B by exponentially amplifying U as above, which has the effect of sharpening the signals in U.
  7. Update the current MEM spectrum F by a forming a weighted average between it and the amplified residual B:
               F[new] = (1 - alpha)*F + alpha*B
    
  8. Repeat steps 1-7 until the chi-square value of the time-domain residual is small enough to be accounted for by the known noise level.

RECONSTRUCTION STRATEGY

In the nmrPipe implementation, the dimensions to be reconstructed by MEM must all be in the same Fourier domain, either as hypercomplex time-domain data, or real-only frequency-domain data. As with linear prediction (LP), the remaining dimensions should usually be in the frequency-domain. Furthermore, phase-correction and baseline adjustment of the spectrum as a whole should be performed prior to using MEM; for this reason, it is usually most convenient to apply MEM to frequency-domain data, rather than time-domain data. So, one useful strategy is to process data completely with the usual methods, inspect the result, and then employ a MEM scheme which will inverse-process and re-process the spectrum in a form most suitable for MEM. A given dimension of an ordinary spectrum might be prepared for MEM with the following steps:

  1. Reconstruct the spectrum's imaginary data.
  2. Undo any phasing previously applied.
  3. Inverse-transform the data into the time-domain.
  4. Undo the zero-filling previously applied.
  5. Remove the window function. (nmrPipe also removes first-point scaling)
  6. Re-apply the first-point scaling.
  7. Zero fill to four times the original size.
  8. Re-transform the data into the frequency-domain.
  9. Re-apply the previous phasing.
These steps are illustrated in the following script fragment:
        | nmrPipe -fn HT -auto                 \
        | nmrPipe -fn PS -inv -hdr             \
        | nmrPipe -fn FT -inv                  \
        | nmrPipe -fn ZF -inv                  \
        | nmrPipe -fn APOD -inv -hdr           \
        | nmrPipe -fn MULT -hdr -xn 1          \
        | nmrPipe -fn ZF -zf 2 -auto           \
        | nmrPipe -fn FT                       \
        | nmrPipe -fn PS -hdr -di              \

Since MEM and MEM schemes involve Hilbert transform and inverse processing, the the following guidelines should be observed:

Zero Filling: in order to reconstruct imaginary data correctly using the Hilbert Transform (HT) during MEM, the original data must have been zero-filled to at least twice its original size. However, additional zero-filling, commonly to four times the original data size, is needed to achieve an increase in resolution by MEM. Furthermore, the transform steps will work much more quickly if the size after zero-fill is a power of two. Note that these zero-filling guidelines will often lead to exceptionally large data sizes.

Acquisition Delay: the Hilbert Transform will only produce ideal reconstruction of imaginary data in two cases: an ordinary HT can be used for data with no acquisition delay (P1=0); a mirror-image HT can be used for data with a half-dwell delay (P1=180). In all other cases, the HT may introduce distortions at the edges of the spectrum, but these may only be serious in cases where the data size is small or where signals of interest lie at the edges of the spectrum.

Window Function: inverse processing will require dividing the data by the original window function. This means that the original window function should be chosen so that none of the window values are close to or equal to zero.

Extracted Regions: the Hilbert Transform, as well as other Fourier manipulations, may not be ideal when applied to extracted regions of a given dimension. This is most likely to be a problem in cases where the data size is small, where signals of interest lie at the edges of the extracted region, or where substantial truncation artifacts extend to the edges of the extracted region.

EXTRACTING AN INDIRECT 2D PLANE

The nmrDraw graphical interface includes commands to interactively select and save 2D planes from the X/Z and Y/Z dimensions of a 3D or 4D spectrum. In typical processing schemes, the directly-acquired dimension will usually be saved as the X-Axis or Y-Axis of the result. In both of these cases, a 2D plane from the indirect dimensions corresponds to an orthogonal plane from nmrDraw's vertical axis, which can be extracted by the following command:

1D Menu/Extract 2D V This nmrDraw command extracts the orthogonal 2D plane which corresponds to the vertical 1D slice currently selected. The 2D plane will be extracted as a file called "ext.dat", which can be treated as a 2D spectrum. The nmrDraw command Previous Data in the File Menu can be used to return to the original 3D data. The Extract 2D command is currently implemented only for all-real data.

Once an appropriate plane is extracted, it can be used to establish and test a MEM scheme and its parameters.

ESTIMATING THE TIME-DOMAIN NOISE LEVEL

Use of MEM requires an estimate of the standard deviation of the noise in the time-domain. Since this may not be easy to measure directly, we can instead use the noise level in the corresponding spectrum to infer the level of noise in the time-domain. The spectral noise level can be estimated interactively with nmrDraw:

  1. Use the "Mouse/2D Zoom" command in nmrDraw to activate the 2D zoom box.
  2. Use the mouse to adjust the size and position of the zoom box to select a small empty region.
  3. Move the mouse pointer outside the graphics area, and click the middle mouse button for "Statistics". A pop-up window will appear, with a histogram of the intensities in the selected region, and a statistical summary. The standard deviation is labeled as "SDev:".

Note that the spectral noise level should be estimated from the data dimensions where MEM will be applied. For example, if 2D MEM will be applied to the CACB/N indirect dimensions of a CBCANH experiment, the noise level in a CACB/N plane should be used, rather than the noise level from an HN/N or HN/CACB plane.

Once the spectral noise level has been estimated, the program showApod can be used to calculate the corresponding time-domain noise level. This program takes into account the details of the window functions and zero filling which were originally used to process the data. For example, if the 2D spectral plane "ext.dat" has a frequency-domain noise level of 9100, the following command could be used:

       showApod -in ext.dat -sigma 9100
to generate a report like the following, which in this case lists the corresponding time-domain noise as roughly 400:
        REMARK Effect of Processing on Noise for ext.dat
        REMARK User-supplied Noise in Processed Data: 9100
        REMARK Noise Before Processing N and CACB: 407.003

        VARS   AXIS LABEL       LW_ADJ  SIGMA_FACTOR
        FORMAT  %8s  %-8s       %+9.5f         %8.5f

         X-Axis     N          0.57981       0.25223
         Y-Axis     CACB       0.56989       0.17732

INTERPRETING THE MEM REPORTS

Setting the report level as -report 2 will produce a summary like the following one, which gives an initial summary of parameters such as memory used, the MEM statistics at each iteration, and a final summary. If the report level is set as -report 1, only the final summary will be listed.

       2D MEM, Two-Channel Mode, 2.65 MBytes
       X-Axis:  1st-Point=Not Adjusted P0=+000.00 P1=+000.00
       Y-Axis:  1st-Point=Not Adjusted P0=-090.00 P1=+180.00
       | 1. Chi/N: 8.75  S: -3.831e+08  La: 0.000  Z: 0.056 |
       | 2. Chi/N: 8.62  S: -3.838e+08  La: 0.151  Z: 0.054 |
       | 3. Chi/N: 7.17  S: -3.953e+08  La: 0.301  Z: 0.059 |
       | 4. Chi/N: 3.16  S: -4.478e+08  La: 0.462  Z: 0.103 |
       | 5. Chi/N: 1.77  S: -4.831e+08  La: 0.704  Z: 0.149 |
       | 6. Chi/N: 1.32  S: -5.029e+08  La: 1.024  Z: 0.181 |
       | 7. Chi/N: 1.13  S: -5.170e+08  La: 1.389  Z: 0.201 |
       | 8. Chi/N: 1.02  S: -5.279e+08  La: 1.779  Z: 0.217 |
       | 9. Chi/N: 0.96  S: -5.362e+08  La: 2.183  Z: 0.227 |
       1. It: 9 C/N: 0.96 S: -5.36e+08 RMS: 4.6e+02 Z: 0.227 CHI2

The meanings of the parameters reported at each iteration are as follows:

Chi/N: Reports the chi2/n value at the current iteration. This should decrease to 1.0 or smaller for a good reconstruction. If this value shows a tendency to increase, try reducing the -alpha parameter or increasing the -eScale parameter by an order of magnitude. If applicable, try using milder deconvolution windows.

If chi2/n converges to a value higher than 1.0, check the accuracy of the noise estimate. If chi2/n value tends to improve too slowly, try increasing the -alpha parameter or decreasing the -eScale parameter.

S: Reports the entropy s at the current iteration. In a good reconstruction. this will be a negative number which increases in magnitude.

A: Reports the value of a, the scaling factor used to create the amplified residual B, which is used in turn to update the current MEM spectrum.

La: Reports the value of lambda at the current iteration. In a good reconstruction, this will be an increasingly large positive number.

Z: Reports another MEM convergence statistic (the cosine of the angle between the entropy gradient and the chi-square gradient), which should ideally increase to 1.0 during the reconstruction. In many practical cases however, this statistic may be much less than 1.0, indicating that the final reconstruction does not have maximal entropy. In these cases, it may be useful to try the scaling mode option -ts 2, which may improve the Z statistics.

PRIMARY OPTIONS

The options listed here, along with the DECONVOLUTION OPTIONS below, are those most likely to be used or adjusted in a typical MEM scheme.

-ndim dimCount [1]
Number of dimensions in each reconstruction; this must be less than or equal to the total number of dimensions in the input data. For example, using a setting of -ndim 2 for a 3D data set means that each 2D XY plane from the 3D data set will be reconstructed separately. Note that memory and computational requirements increase rapidly with increasing dimension count.

-sigma noise [1.0]
This important parameter specifies the estimated standard deviation of the noise in the time-domain. If this value is not set well, a poor reconstruction will probably result. The time-domain noise level is usually not measured directly, but rather it is estimated on the basis of the noise level measured in the frequency-domain. See the section on Noise Estimation above for more. If a -sigma value of zero is given, an automated estimate will be used, but this is not recommended. In most cases, it seems better to underestimate the noise level, in order to prevent MEM iterations from terminating before the smallest signals have been reconstructed. If the MEM result seems to be missing signals that are observed in the original data, try reducing the value of the -sigma parameter.

-alpha alpha [1e-3]
In each iteration of the MEM algorithm, an "ideal" correction is calculated for updating the current MEM spectrum. In practice however, only a fraction of this correction is added to insure stable convergence. This option specifies the parameter alpha, the fraction of the ideal correction to add to the current MEM spectrum. Values for alpha commonly range from 1e-1 to 1e-4, with larger values leading to faster but potentially unstable reconstructions. Many divergence problems can be fixed by reducing the alpha parameter.

-x0 thresh [5.0]
In the MEM algorithm, the scaled frequency-domain residual (2*lambda*U) at a given iteration is amplified exponentially to generate the correction to add to the current MEM spectrum. To avoid excessive amplification, scaled signals above the threshold given by -x0 will be amplified linearly rather than exponentially. Increasing this parameter may improve resolution or convergence speed, at the expense of stability.

-eScale g [1.0]
Specifies the scale factor g used for computing the entropy s from MEM spectrum F, s = -Sum( F ln(F/g) ). Smaller values will lead to faster but potentially unstable reconstructions.

-ts tsMode [0]
This option specifies when to perform additional scaling of the current MEM spectrum in an attempt to more effectively minimize the time-domain residual. The possible modes are:
 0 No scaling is used (Default).
 1 Data is rescaled on divergence.
 2 Data is rescaled at every iteration.

Option 2 can often improve the convergence of the Z statistic, possibly at the expense of increased iteration count.

-report rLevel [0]
Specifies the report level for statistics displayed during the calculation:
 0 No report (Default).
 1 Report summary after each final iteration.
 2 Report statistics after each iteration
 3 Save current data F and U at each iteration.

Option 3 is intended primarily for testing and development purposes.

DECONVOLUTION OPTIONS

These options specify the name and parameters of the window functions used on the dimensions of the MEM time-domain data (R in the description above) before it is compared with the original time-domain data. These windows will usually be line broadening functions whose broadening widths are less than the intrinsic linewidths in the original data. The purpose of these windows is to adjust the (sharp) MEM spectrum so that it matches the original data when broadened. In other words, resolution enhancement is achieved via a convolution of the MEM data (which is a stable procedure), rather than by direct deconvolution of the original data (which is unstable).

A given window is specified as a function name option, plus one, two, or three generic window parameters, called Q1, Q2, and Q3. The dimension that the window and parameters apply to is specified by including a prefix "x", "y", "z", or "a" in the option name, such as -xconv, -yconv, etc. If the prefix is omitted, such as in -conv, the option applies to the X-Axis of the data.

The specified function name can be any nmrPipe window function, including SP, EM, or GM. Depending on which window function is selected, the corresponding options -cQ1, -cQ2, etc. specify its Q parameters; the meaning of the parameters depends on the choice of window function. The specific interpretation of the Q1/Q2/Q3 parameters is given in the manual page about a given window function, and also in the nmrPipe on-line help text, for example by using:

        nmrPipe -help -fn SP
The complete list of deconvolution options follows:

-conv  cFn  (for X-Axis Only)
-xconv xcFn (for X-Axis) -yconv ycFn (for Y-Axis) -zconv zcFn (for Z-Axis)
Specifies the window function name for the given axis of the MEM time-domain data; if no option is specified, no window will be applied.

-cQ1 cQ1 -cQ2 cQ2 -cQ3 cQ3
-xcQ1 xcQ1 -xcQ2 xcQ2 -xcQ3 xcQ3
Specifies the window parameters Q1, Q2, and Q3 for the window applied to the X-Axis of the MEM time-domain data. These options are only used if options -conv or -xconv have been used as well.

-ycQ1 ycQ1 -ycQ2 ycQ2 -ycQ3 ycQ3
Specifies the window parameters Q1, Q2, and Q3 for the window applied to the Y-Axis of the MEM time-domain data. These options are only used if option -yconv has been used as well.

-zcQ1 zcQ1 -zcQ2 zcQ2 -zcQ3 zcQ3
Specifies the window parameters Q1, Q2, and Q3 for the window applied to the Z-Axis of the MEM time-domain data. These options are only used if option -zconv has been used as well.

CONVERGENCE TEST OPTIONS

-cLim chiLimit [1.0]
Probably the most important convergence criterion for MEM is that the residual between the original time-domain data and the MEM reconstruction is small enough to be accounted for by the known noise level. This is judged by the chi-square value of the time-domain residual. Ideally, for n total points in the data (chi-square)/n should be one or smaller. The option -cLim specifies the value which (chi-square)/n should fall below in a properly converged result. Iterations will terminate automatically if this criterion is reached.

-zLim zLimit [1.0]
One convergence criterion for MEM is the balance between corrections which maximize the entropy of the current MEM spectrum compared to corrections which minimize the chi-square value of the residual. This is measured by the MEM Z value, which gives the cosine of the angle between the entropy gradient and the chi-square gradient. Iterations will terminate automatically if the Z value is greater than the given zLimit parameter. Note that the Z value should increase to 1.0 in an ideal reconstruction, but values as small as 0.1 often seem to give acceptable results.

-itmax maxIter [50]
Specifies the maximum iteration count allowed for any given reconstruction. Values in the range of 10 to 100 are usual.

-itmin minIter [3]
Specifies the minimum iteration count allowed for any given reconstruction, even if convergence tests have already been satisfied.

-div maxDiverge [3]
Specifies the maximum number of successive iterations where the chi-square value is permitted to get worse (larger) before iterations terminate automatically. Note that divergence is often due to an -alpha parameter which is too large, a bad noise estimate for -sigma, or a problem in the pre-processing scheme.

RECONSTRUCTION MODE OPTIONS

-neg This flag enables reconstruction of combinations of both positive and negative signals (Two Channel MEM) It is the default mode.

-pos This flag limits MEM to reconstruction of only positive signals. If the data contains substantial negative signals, this mode may cause a poor result and bad convergence statistics.

-zero The MEM algorithm used here tends to introduce a zero-order offset into the reconstruction. This option will cause the final result to be corrected by zero-order baseline estimation. (Default).

-nozero This flag will suppress the zero-order offset correction which is applied to the final result.

-freq This flag will produce the final MEM result as frequency-domain data. (Default).

-time This flag will produce the final MEM result as time-domain data.

-zf zfCount [2] (Also -xzf -yzf etc).
Specifies the zero fill count for a given dimension. The zero fill count defines the number of times to double the data size by padding with zeros. Because the MEM algorithm used here includes a Hilbert transform step, each dimension in the reconstruction should be zero-filled at least once for MEM to work correctly, and optimally twice (i.e. zfCount of 2, to increase size by a factor of 4) to make resolution improvement possible. Note however that it is recommended that zero-filling is performed as a pre-processing step prior to use of MEM. See the section on Reconstruction Schemes above for more.

OTHER MEM PARAMETER OPTIONS

-clip cMode [0]
The option activates clipping of the original spectrum, to reduce the most extreme intensities to a given level. It can be used to stabilize reconstruction of data with very large peaks which are not of direct interest (e.g. solvent signal). The possible modes are:

 0 No clipping (Default).
 1 Clip positive signals only.
 2 Clip negative signals only.
 3 Clip both positive and negative signals.

-min cMin [-1e17]
Specifies the minimum value for negative signal clipping; this only applies for -clip modes 2 and 3 above.

-max cMax [+1e17]
Specifies the minimum value for negative signal clipping; this only applies for -clip modes 1 and 3 above.

-cScale x [1e2]
Specifies a scaling parameter x; the original time-domain data will be multiplied by x/sigma at the start of the MEM calculation, and the final MEM result will be multiplied by sigma/x to remove the effect of the initial scaling. The intention of this scaling is to keep spectral intensities within a reasonable range during the reconstruction.

-lamb lambda [1e-5]
Specifies the initial value for the La Grange multiplier lambda; this value is also used as the minimum allowable value for lambda at any given iteration.

-step beta [0.02]
At each iteration, an ideal value of the the La Grange multiplier lambda is computed. However this "ideal" calculation can lead to lambda values which change too quickly from one iteration to the next, leading to instability and divergence. The beta parameter specified here is used to limit the change in the lambda value according to:

             lambda[next] = (1 - beta)*lambda[prev] + beta*lambda[ideal]
Larger values of beta will lead to faster but potentially unstable reconstructions.

-tScale w [0.5]
Specifies a scale factor w for enabling exponential amplification, as a value in the range of 0 to 1. If this value is non-zero, only data points larger than w*Max( U ) will be amplified, where U is the frequency-domain residual between the current MEM spectrum and the original data. The intent of this parameter is to avoid amplification of truncation artifacts. Note however that this improvement may be at the expense of lineshape distortion.

-tLimit tLimit [1.5]
Specifies a chi2 limit for use of the amplification threshold specified by -tScale w. When the chi2 value goes below this limit, the -tScale option will no longer suppress amplification.

EXAMPLES

1D MEM applied to a 2D Spectrum: MEM may sometimes be useful in reducing truncation artifacts or improving resolution in 2D spectra when applied separately to each 1D t1 vector. Note however that differences in convergence between adjacent 1D vectors may distort the overall 2D peak shapes. The following scheme was applied to a 2D proton-carbon HSQC spectrum:

      #!/bin/csh

      nmrPipe -in test.fid \
      | nmrPipe  -fn SP -off 0.5 -end 0.95 -pow 1 -c 0.5 \
      | nmrPipe  -fn ZF -auto                            \
      | nmrPipe  -fn FT -auto                            \
      | nmrPipe  -fn PS -p0 11 -p1 -22 -di               \
      | nmrPipe  -fn EXT -x1 6ppm -xn 2ppm -sw           \
      | nmrPipe  -fn TP                                  \
      | nmrPipe  -fn MULT -c 0.5 -xn 1                   \
      | nmrPipe  -fn ZF -zf 2 -auto                      \
      | nmrPipe  -fn FT -auto                            \
      | nmrPipe  -fn PS -p0 0.0 -p1 0.0 -di              \
      | nmrPipe  -fn MEM -sigma 3000 -report 2           \
      | nmrPipe  -fn TP                                  \
         -ov -out mem1d.ft2

2D MEM applied to an Extracted 2D Plane: the following macro gives a scheme for re-processing a 2D plane extracted from a 3D spectrum. The scheme expects that the dimensions of the 2D plane where originally processed in such a way as to permit inverse processing.

      #!/bin/csh

      nmrPipe -in ext.dat \
      | nmrPipe -fn HT -auto                           \
      | nmrPipe -fn PS -inv -hdr                       \
      | nmrPipe -fn FT -inv                            \
      | nmrPipe -fn ZF -inv                            \
      | nmrPipe -fn SP -inv -hdr                       \
      | nmrPipe -fn MULT -hdr -xn 1                    \
      | nmrPipe -fn ZF -zf 2 -auto                     \
      | nmrPipe -fn FT                                 \
      | nmrPipe -fn PS -hdr -di                        \
      | nmrPipe -fn TP                                 \
      | nmrPipe -fn HT -auto                           \
      | nmrPipe -fn PS -inv -hdr                       \
      | nmrPipe -fn FT -inv                            \
      | nmrPipe -fn ZF -inv                            \
      | nmrPipe -fn SP -inv -hdr                       \
      | nmrPipe -fn MULT -hdr -xn 1                    \
      | nmrPipe -fn ZF -zf 2 -auto                     \
      | nmrPipe -fn FT                                 \
      | nmrPipe -fn PS -hdr -di                        \
      | nmrPipe -fn TP                                 \
      | nmrPipe -fn MEM -ndim 2 -sigma 460 -report 2   \
         -out ext.mem.ft2 -ov

2D MEM applied to a 3D Spectrum: the following macro gives a scheme for re-processing a 3D spectrum with 2D MEM in both indirect dimensions. This macro expects that the original spectrum is stored in transposed order, with the directly-detected dimension in the Y-Axis. It also expects that the indirect dimensions where originally processed in such a way as to permit inverse processing.

The scheme has been built to present the data to MEM in a specific axis order. For instance, if we assume the input spectrum has a data order X=H Y=HN Z=N, then this scheme temporarily re-arranges the data for MEM such that X=N Y=H Z=HN.

In the macro, the indirect planes are reprocessed so that they are presented for MEM without window functions, but with first-point scaling in the time-domain, phasing, and extensive zero fill. MEM deconvolution is used to apply exponential line sharpening (EM) by 10hz in the X-Axis, and 15hz in the Y-Axis:

       #!/bin/csh

       xyz2pipe -in ft/test%03d.ft3 -z                  \
       | nmrPipe -fn HT -auto                           \
       | nmrPipe -fn PS -inv -hdr                       \
       | nmrPipe -fn FT -inv                            \
       | nmrPipe -fn ZF -inv                            \
       | nmrPipe -fn SP -inv -hdr                       \
       | nmrPipe -fn MULT -c 0.5 -xn 1                  \
       | nmrPipe -fn ZF -zf 2 -auto                     \
       | nmrPipe -fn FT                                 \
       | nmrPipe -fn PS -hdr -di                        \
       | nmrPipe -fn TP                                 \
       | nmrPipe -fn HT -auto                           \
       | nmrPipe -fn PS -inv -hdr                       \
       | nmrPipe -fn FT -inv                            \
       | nmrPipe -fn ZF -inv                            \
       | nmrPipe -fn SP -inv -hdr                       \
       | nmrPipe -fn MULT -c 0.5 -xn 1                  \
       | nmrPipe -fn ZF -zf 2 -auto                     \
       | nmrPipe -fn FT                                 \
       | nmrPipe -fn PS -hdr -di                        \
       | nmrPipe -fn TP                                 \
       | nmrPipe -fn MEM -ndim 2 -report 2              \
           -sigma 1520 -alpha 0.0001 -eScale 0.1        \
           -xconv EM -xcQ1 10 -yconv EM -ycQ1 15        \
       | pipe2xyz -out ft/mem%03d.ft3 -z

Deconvolution of a fixed coupling: in this example, a J = 88Hz coupling is deconvolved from the indirect dimension of 2D time-domain data. The coupling profile is specified via parameters for the JMOD window function. In the NMRPipe implementation of MEM, this kind of deconvolution will only work effectively on data with limited dynamic range. The calculation is stabilized by using a higher than normal threshold for computing an iterate (-tScale 0.8). Furthermore, in this case, the iterations are stopped at an earlier than usual stage (-cLim 3.0) to avoid artifacts from the deconvolution, which take the form of satellite peaks at positions +/- J from the deconvolved peaks.

In this example, the data are temporarily transformed in the direct dimension in order to extract a PPM range of interest.

nmrPipe -in test.fid \
| nmrPipe -fn ZF -auto                             \
| nmrPipe -fn FT                                   \
| nmrPipe -fn EXT -x1 12ppm -xn 6ppm -sw           \
| nmrPipe -fn PS -p0 256 -p1 0                     \
| nmrPipe -fn FT -inv \
| nmrPipe -fn ZF -inv \
| nmrPipe -fn MEM -sigma 750 -report 2 -ndim 2 -cLim 3 -tScale 0.8 -tLimit 0.0 \
          -xconv EM   -xcQ1 10.0 \
          -yconv JMOD -ycQ1  0.5 -ycQ2 88.0 -ycQ3 5.0 -yzf 3 \
          -out deco.ft2 -ov 
HEADER VALUES

MEM uses NDFTFLAG to determine the transform state of the current data. It also uses NDP1 and NDX1/NDXN to decide whether an ordinary Hilbert Transform should be used to reconstruct imaginary data, or a special version for data with half-point delay (NDP1 = 180, no sub-region extracted).

MEM updates the data sizes NDSIZE, NDAPOD, and NDTDSIZE. The result from MEM is real-only frequency domain data, so that MEM also updates NDFTFLAG, NDQUADFLAG, and FDQUADFLAG, as well as the parameters for PPM calibration (NDORIG, NDCENTER).