*! ralslm 1.0.1 16may2026 Dr Merwan Roudane *! RALS-LM unit-root test of Meng, Im, Lee & Tieslau (2014) * Festschrift in Honor of Peter Schmidt (pp. 343-357), Springer. * Mirrors the GAUSS code rals_lm.src by Saban Nazlioglu. *------------------------------------------------------------------------------ program define ralslm, rclass version 14.0 capture mata: __rals_loaded() if _rc qui _rals_mata syntax varname(ts) [if] [in], [ /// MAXLags(integer 8) /// IC(string) /// Graph /// noHEADer ] marksample touse qui count if `touse' if r(N) < 25 { di as error "sample too small to compute RALS-LM (N=" r(N) ")" exit 2001 } capture tsset if _rc { di as error "data must be {bf:tsset} before running ralslm" exit 459 } local timevar `r(timevar)' if "`ic'"=="" local ic "tstat" local ic = lower("`ic'") if !inlist("`ic'", "aic", "bic", "tstat") { di as error "ic() must be aic, bic, or tstat" exit 198 } tempname results cv mata: _rals_lm_engine("`varlist'", "`touse'", `maxlags', "`ic'", "`results'", "`cv'") local tauLM = `results'[1,1] local tauRALS = `results'[1,2] local rho2 = `results'[1,3] local T = `results'[1,4] local lag = `results'[1,5] if "`header'"=="" { _rals_print_header "RALS-LM" "Meng, Im, Lee & Tieslau (2014)" "`varlist'" `T' `lag' "`ic'" 2 } _rals_print_lm_table `tauLM' `tauRALS' `rho2' `cv' return scalar tauLM = `tauLM' return scalar tauRALS = `tauRALS' return scalar rho2 = `rho2' return scalar T = `T' return scalar lag = `lag' return scalar cv01_LM = `cv'[1,1] return scalar cv05_LM = `cv'[1,2] return scalar cv10_LM = `cv'[1,3] return scalar cv01 = `cv'[2,1] return scalar cv05 = `cv'[2,2] return scalar cv10 = `cv'[2,3] return local ic = "`ic'" return local test = "RALS-LM" return local cmdline = "ralslm `varlist', maxlags(`maxlags') ic(`ic')" if "`graph'" != "" { _rals_graph_ur "`varlist'" "`timevar'" "`touse'" "RALS-LM" `tauRALS' `cv'[2,2] } end program define _rals_print_lm_table args tauLM tauRALS rho2 cv _rals_print_two_row "LM (stage 1)" "RALS-LM (stage 2)" /// `tauLM' `tauRALS' `rho2' /// `cv'[1,1] `cv'[1,2] `cv'[1,3] `cv'[2,1] `cv'[2,2] `cv'[2,3] local concl_lm = cond(`tauLM' < `cv'[1,2], "{bf:reject} H0 (LM, 5%)", /// "fail to reject H0 (LM, 5%)") local concl_rals = cond(`tauRALS' < `cv'[2,2], "{bf:reject} H0 (RALS-LM, 5%)", /// "fail to reject H0 (RALS-LM, 5%)") di as text " Decision (5%)" di as text " {c -}{c -} LM : " as result "`concl_lm'" di as text " {c -}{c -} RALS-LM : " as result "`concl_rals'" di as text "{hline 80}" di as text " Stage-1 CVs : Schmidt-Phillips (1992) trend-LM table." di as text " Stage-2 CVs : Meng, Im, Lee & Tieslau (2014) rho^2 interpolation." di as text "" end version 14.0 mata: mata set matastrict off void _rals_lm_engine(string scalar yname, string scalar touse, real scalar pmax, string scalar icmode, string scalar rname, string scalar cvname) { real colvector y, dy, dz, b0, ylm, ly, dep, ds, depf, es, ea real scalar T, p, psi, sig2, sig2A, tauLM, tauRALS, m2, m3, rho2 real matrix lmat, Xs, Xa, XXi, XXi2, w, crt, cvr real colvector bs, ses, ba, sea real scalar popt, bestic, j st_view(y, ., yname, touse) T = rows(y) dy = y[2::T] - y[1::T-1] // Schmidt-Phillips style detrending -------------------------------------- // dy = b0 + e ; b0 is the OLS slope on a constant, i.e. the mean of dy. real scalar b0s b0s = sum(dy)/rows(dy) psi = y[1] - b0s // detrended series real colvector ylm0 ylm0 = y :- psi :- (1::T)*b0s ylm = ylm0[2::T] // drop first as in src real scalar Tlm Tlm = rows(ylm) ly = ylm[1::Tlm-1] // y_{t-1} ds = ylm[2::Tlm] - ylm[1::Tlm-1] // Δ of detrended // dy needs same trimming real colvector dyt dyt = dy[2::rows(dy)] // align with ds // build lag matrix of ds lmat = J(rows(ds), pmax, 0) for (j=1; j<=pmax; j++) lmat[j+1::rows(ds), j] = ds[1::rows(ds)-j] // full design at pmax: [ly~const~ldy], match rals_lm.src real matrix Xfull Xfull = ly, J(rows(ly),1,1), lmat real colvector depfull depfull = dyt if (pmax>0) { depfull = depfull[pmax+1::rows(dyt)] Xfull = Xfull[pmax+1::rows(dyt), .] } // lag selection ---------------------------------------------------------- real colvector icvec, tvec icvec = J(pmax+1,1,.); tvec = J(pmax+1,1,.) for (p=0; p<=pmax; p++) { real matrix Xp; real matrix XXp; real colvector bpv,epv,sep Xp = Xfull[., 1::(2+p)] XXp = invsym(quadcross(Xp,Xp)) bpv = XXp*quadcross(Xp,depfull) epv = depfull - Xp*bpv if (icmode=="bic") { icvec[p+1] = ln(epv'epv/rows(Xp))+cols(Xp)*ln(rows(Xp))/rows(Xp) } else { icvec[p+1] = ln(epv'epv/rows(Xp))+2*cols(Xp)/rows(Xp) } if (p>0) { sep = sqrt(diagonal((epv'epv/(rows(Xp)-cols(Xp)))*XXp)) tvec[p+1] = abs(bpv[cols(Xp),1]/sep[cols(Xp),1]) } } if (icmode=="tstat") { popt = 0 for (p=pmax; p>=1; p--) { if (tvec[p+1]>=1.645) { popt = p break } } } else { popt = 0 bestic = icvec[1] for (p=1; p<=pmax; p++) { if (icvec[p+1] < bestic) { bestic = icvec[p+1] popt = p } } } // stage 1: LM stat at chosen lag ----------------------------------------- Xs = Xfull[., 1::(2+popt)] XXi = invsym(quadcross(Xs,Xs)) bs = XXi*quadcross(Xs,depfull) es = depfull - Xs*bs sig2 = (es'es)/(rows(Xs)-cols(Xs)) ses = sqrt(diagonal(sig2*XXi)) tauLM = bs[1]/ses[1] // RALS augmentation ------------------------------------------------------ m2 = sum(es:^2)/rows(es) m3 = sum(es:^3)/rows(es) w = (es:^2 :- m2), (es:^3 :- m3 :- 3*m2*es) Xa = Xs, w XXi2 = invsym(quadcross(Xa,Xa)) ba = XXi2*quadcross(Xa,depfull) ea = depfull - Xa*ba sig2A = (ea'ea)/(rows(Xa)-cols(Xa)) sea = sqrt(diagonal(sig2A*XXi2)) tauRALS = ba[1]/sea[1] rho2 = sig2A/sig2 if (rho2<0) rho2=0 if (rho2>1) rho2=1 crt = __rals_cv_lm() cvr = __rals_interp(crt, rho2) real rowvector cvLM cvLM = __rals_cv_lmstat(T-popt) real matrix CV2 CV2 = cvLM \ cvr st_matrix(rname, (tauLM, tauRALS, rho2, T, popt)) st_matrix(cvname, CV2) } end