mtmML96(dat,tbw=3,ntap=NULL,padfac=5,demean=T,detrend=F,medsmooth=0.2,
opt=1,linLog=2,siglevel=0.9,output=0,CLpwr=T,xmin=0,xmax=Nyq,
sigID=F,pl=1,genplot=T,verbose=T)
To help address this issue, an alternative median smoothing approach is applied that implements Tukey's robust end-point rule and symmetrical medians (see the function runmed for details). Numerical experiments indicate that this approach produces an approximately uniform false positive rate across the spectrum. It should be noted that the false positive rates are still inflated with this method, but they are substantially reduced compared to the original ML96 approach. For example, simulations using rho=0.9 (using identical parameters to those in Meyers, 2012) yield median false positive rates of 1.7%, 7.3% and 13.4%, for the 99%, 95% and 90% confidence levels (respectively). This compares with 4.7%, 11.4% and 17.8% using the original approach (see Table 2 of Meyers, 2012).
NOTE: If the (fast) Brent or Gauss-Newton methods fail, use the (slow) grid search approach.
This version of the ML96 algorithm was first implemented in Patterson et al. (2014).
Meyers, S.R., 2012, Seeing red in cyclic stratigraphy: Spectral noise estimation for astrochronology, Paleoceanography, 27, PA3228.
Patterson, M.O., McKay, R., Naish, T., Escutia, C., Jimenez-Espejo, F.J., Raymo, M.E., Meyers, S.R., Tauxe, L., Brinkhuis, H., and IODP Expedition 318 Scientists, 2014, Response of the East Antarctic Ice Sheet to orbital forcing during the Pliocene and Early Pleistocene, Nature Geoscience, v. 7, p. 841-847.
Thomson, D. J., 1982, Spectrum estimation and harmonic analysis, Proc. IEEE, 70, 1055-1096, doi:10.1109/PROC.1982.12433.
http://www.meteo.psu.edu/holocene/public_html/Mann/tools/MTM-RED/
Tukey, J.W., 1977, Exploratory Data Analysis, Addison.
runmed
, spec.mtm
, mtmAR
, lowspec
, and periodogram
# generate example series with periods of 400 ka, 100 ka, 40 ka and 20 ka
ex = cycles(freqs=c(1/400,1/100,1/40,1/20),start=1,end=1000,dt=5)
# add AR1 noise
noise = ar1(npts=200,dt=5,sd=0.5)
ex[2] = ex[2] + noise[2]
# run ML96 analysis
pl(1, title="mtmML96")
mtmML96(ex)
# compare to analysis with conventional AR1 noise test
pl(1,title="mtm")
mtm(ex)
# compare to analysis with LOWSPEC
pl(1, title="lowspec")
lowspec(ex)
# compare to amplitudes from eha
pl(1,title="eha")
eha(ex,tbw=3,win=1000,pad=1000)
Run the code above in your browser using DataLab