:Python分段建模時間序列突變)
簡介門限自回歸TAR模型是在傳統(tǒng)線性AR模型基礎(chǔ)上引入閾值機制的非線性時間序列方法適用于刻畫不同狀態(tài)下回歸關(guān)系發(fā)生結(jié)構(gòu)性變化的數(shù)據(jù)。這份MATLAB代碼資源面向時間序列分析方向的學(xué)習(xí)者與研究人員圍繞TAR模型構(gòu)建與模型選擇提供了一套可運行的完整示例。壓縮包共6個文件含3個m腳本和3個txt文件體積僅70KB腳本覆蓋數(shù)據(jù)預(yù)處理、閾值識別、分段AR擬合、最大似然參數(shù)估計及LR似然比圖繪制等核心環(huán)節(jié)txt文件提供示例數(shù)據(jù)與使用說明便于對照運行。通過調(diào)試和運行這套代碼讀者可掌握TAR模型從建模到診斷的完整流程并學(xué)會利用LR圖判斷最佳門限數(shù)量提升對非線性時間序列的建模能力。目前已有356人學(xué)習(xí)下載適合具備基礎(chǔ)回歸知識、希望實踐門限模型的MATLAB用戶。1. 門限自回歸時間序列回歸里的“變臉”問題終于有解了做時間序列預(yù)測的人遲早會撞上一堵墻模型在某個時間點之前擬合得漂亮之后卻集體失效。不是數(shù)據(jù)變臟了而是數(shù)據(jù)背后的“機制”變了——宏觀數(shù)據(jù)遇到政策拐點風(fēng)速序列遇到季節(jié)轉(zhuǎn)換銷售數(shù)據(jù)遇到大促前后。線性 AR 模型假設(shè)過去的關(guān)系永遠(yuǎn)不變但現(xiàn)實里很多序列存在“門限效應(yīng)”當(dāng)某個驅(qū)動變量跨過閾值序列的行為就會切換成另一套邏輯。門限自回歸Threshold Autoregressive Model簡稱 TAR就是專門為這種“變關(guān)系”設(shè)計的模型。它不強行用一條曲線去擬合所有時期而是把序列按門限變量的取值切成若干區(qū)間每個區(qū)間內(nèi)分別擬合自回歸方程。jasa_03m 這類月度數(shù)據(jù)正是 TAR 模型最典型的用武之地月度觀測里天然藏著季節(jié)、周期和突變點用門限模型能比常規(guī) ARIMA 多抓住一層結(jié)構(gòu)。本文不繞彎子直接講清楚門限自回歸的建模流程、Python 實現(xiàn)和參數(shù)設(shè)置給你一份能照著復(fù)現(xiàn)的實戰(zhàn)路徑。2. 門限自回歸的核心邏輯為什么要“分段”做回歸2.1 線性 AR 的假設(shè)盲區(qū)關(guān)系不是恒定的先看一個最基礎(chǔ)的問題門限自回歸到底在解決什么。經(jīng)典的 AR(p) 模型寫成 y_t c φ?y_{t-1} … φ_p y_{t-p} ε_t它隱含的假設(shè)是從第一期到最后一期φ 這些系數(shù)是不變的。換句話說無論序列處于高位、低位還是中間狀態(tài)過去值對當(dāng)前值的影響力度都一樣。這個假設(shè)在平穩(wěn)序列上還能應(yīng)付一旦序列存在結(jié)構(gòu)性變化——比如某個月出臺新政策、市場需求突然換擋——單一組的 φ 就完全不夠用了。門限自回歸的破局思路很直接把序列按某個變量的數(shù)值劃分成不同的“體制”regime每個體制內(nèi)各跑各的自回歸。最經(jīng)典的二體制 TAR 模型寫成y_t (φ?? φ??y_{t-1} … φ??y_{t-p}) · I(z_t ≤ c) (φ?? φ??y_{t-1} … φ??y_{t-p}) · I(z_t c) ε_t其中 z_t 是門限變量c 是門限值I(·) 是指示函數(shù)。公式看著復(fù)雜意思是當(dāng) z_t 小于等于 c 時序列服從第一套自回歸方程當(dāng) z_t 大于 c 時切換到第二套方程。z_t 可以是滯后值 y_{t-d}這叫 SETAR自激勵門限自回歸也可以是外生變量這叫 TAR。這個分段機制讓模型天然具備“識別變結(jié)構(gòu)”的能力。2.2 TAR 與 SETAR 怎么選門限變量是核心分水嶺建模前必須先想清楚你的門限變量是什么。這個問題決定了模型的類別。如果門限變量是序列自身的滯后值比如 z_t y_{t-1}模型就叫 SETAR如果門限變量是另一個外生變量比如利率、匯率、氣溫模型就是 TAR。兩者的共性和差異可以用一個判斷框架說清楚。維度SETARTAR門限變量序列自身滯后值 y_{t-d}外生變量 z_t適用場景序列有自我強化的周期或波動聚集序列受外部變量驅(qū)動切換機制建模復(fù)雜度較低數(shù)據(jù)本身夠用需要配套的門限變量數(shù)據(jù)典型應(yīng)用太陽黑子、氣溫、股價波動銷售額受政策利率影響、風(fēng)速受氣壓影響我個人的習(xí)慣是先跑一遍 SETAR 做摸底因為不需要額外找數(shù)據(jù)如果檢驗出明顯的門限效應(yīng)但不穩(wěn)定再考慮換外生門限變量。見過不少項目是數(shù)據(jù)明明有突變點卻硬套 ARIMA 硬扛改造成門限模型后誤差直接降了一截。選擇 TAR 還是 SETAR本質(zhì)是在問你更相信“過去的狀態(tài)”還是“外部的力量”能解釋當(dāng)前的行為切換。2.3 門限效應(yīng)的檢驗不要上來就分段這里有個新手常犯的錯誤看到序列圖上有幾個折點就急著定門限值。門限模型的第一個正式步驟應(yīng)該是假說檢驗——先驗證“是否存在門限效應(yīng)”再估計門限值。統(tǒng)計上一般用 Chan (1993) 的超檢驗過程先按門限變量的取值排序?qū)γ恳粋€可能作為門限的候選點把樣本分成兩段分別估計兩段的自回歸模型算殘差平方和RSS取 RSS 最小的那個候選點作為門限估計值。這一步在 R 的 tsDyn 包和 Python 里都有成熟實現(xiàn)。但這還不夠因為哪怕沒有真正的門限效應(yīng)你也能找到“使 RSS 最小”的點只是它沒有統(tǒng)計顯著意義。標(biāo)準(zhǔn)的做法是用自舉bootstrap生成零分布在原假設(shè)“無門限效應(yīng)”下模擬一大批序列計算每個模擬序列的超檢驗統(tǒng)計量然后看真實統(tǒng)計量落在分布的多極端位置算出 p 值。如果 p 值小于 0.05才算有統(tǒng)計依據(jù)做分段回歸。3. 用 Python 實現(xiàn)門限自回歸從數(shù)據(jù)檢驗到參數(shù)估計3.1 數(shù)據(jù)準(zhǔn)備與平穩(wěn)性預(yù)檢驗沒有一個可靠的過程就沒有一個可靠的模型任何門限模型的起點都是數(shù)據(jù)清洗和平穩(wěn)性檢驗。不像線性回歸可以容忍一些粗糙的預(yù)處理門限模型對數(shù)據(jù)的“分層”非常敏感如果序列里混著異常值或趨勢項門限點會被強行拉到一個錯誤的位置。先做 ADF 檢驗確認(rèn)平穩(wěn)性不平穩(wěn)就差分這步?jīng)]有繞過空間。import numpy as np import pandas as pd from statsmodels.tsa.stattools import adfuller # 加載 jasa_03m 月度數(shù)據(jù)假設(shè)為單列時間序列 series pd.read_csv(jasa_03m.csv, index_col0, parse_datesTrue).iloc[:, 0] # ADF 檢驗p 值小于 0.05 視為平穩(wěn) adf_result adfuller(series.dropna()) print(fADF 統(tǒng)計量: {adf_result[0]:.4f}, p 值: {adf_result[1]:.4f}) # 若 p 值偏大做一階差分后重新檢驗 if adf_result[1] 0.05: series_diff series.diff().dropna() print(原序列非平穩(wěn)使用一階差分序列) else: series_diff series代碼的邏輯很直觀先用adfuller做 ADF 單位根檢驗p 值小說明沒有單位根、序列平穩(wěn)p 值大就做一階差分。這里有個參數(shù)注意點adfuller默認(rèn)的回歸項含常數(shù)項c如果序列均值明顯非零這個設(shè)置是對的如果序列圍繞 0 波動可以傳入regressionn去掉常數(shù)項能提高檢驗功效。還有一點容易忽略——月度數(shù)據(jù)應(yīng)該把autolagAIC設(shè)為默認(rèn)的自動滯后選擇它會按信息準(zhǔn)則挑最優(yōu)滯后階數(shù)比固定滯后更穩(wěn)。3.2 門限候選點的網(wǎng)格搜索把連續(xù)閾值變成可計算的離散候選門限值本身是連續(xù)的但實際計算時只能從觀測值里選。常見做法是把門限變量 z_t 從小到大排序然后按一定的分位數(shù)區(qū)間比如 15% 到 85%作為候選范圍去掉兩端太少的樣本保證每個體制內(nèi)的觀測數(shù)足夠做回歸。接下來對每個候選點把全樣本切兩半分別估計 AR 模型累加兩個子模型的殘差平方和。import itertools from statsmodels.tsa.ar_model import AutoReg def fit_tar_rss(series, delay, threshold_candidates, ar_order): 遍歷候選門限點返回 RSS 最小的門限值及對應(yīng)模型信息 best {rss: np.inf, threshold: None, models: None} for th in threshold_candidates: # 按門限變量這里是 y_{t-delay}分割樣本 mask_low series.shift(delay) th mask_high series.shift(delay) th s_low, s_high series.loc[mask_low].dropna(), series.loc[mask_high].dropna() # 每個體制內(nèi)要求最小樣本量否則跳過 if len(s_low) 20 or len(s_high) 20: continue model_low AutoReg(s_low, lagsar_order).fit() model_high AutoReg(s_high, lagsar_order).fit() rss_total model_low.ssr model_high.ssr if rss_total best[rss]: best[rss] rss_total best[threshold] th best[models] (model_low, model_high, mask_low, mask_high) return best # 以滯后 1 期為門限變量候選門限設(shè)為序列 15%-85% 分位區(qū)間內(nèi)的觀測值 series_clean series_diff.dropna() threshold_candidates series_clean.shift(1).dropna().quantile([0.15, 0.25, 0.5, 0.75, 0.85]).values result fit_tar_rss(series_clean, delay1, threshold_candidatesthreshold_candidates, ar_order2) print(f最優(yōu)門限值: {result[threshold]:.4f}, 最小 RSS: {result[rss]:.4f})這段代碼的核心是fit_tar_rss函數(shù)里的雙重循環(huán)邏輯外層遍歷門限候選點內(nèi)層對分段后的兩個子樣本分別擬合AutoReg模型。model_low.ssr和model_high.ssr分別是低體制和高體制模型的殘差平方和加起來就是總 RSS。延遲參數(shù)delay決定門限變量用的是哪一期滯后——這個參數(shù)直接影響模型行為后面會專門講怎么選。我在項目里的經(jīng)驗是候選門限別用手工指定用我這段的分位數(shù)切片方式更客觀5 個分位點做粗篩找到最優(yōu)區(qū)間后在區(qū)間內(nèi)細(xì)化網(wǎng)格重跑一遍能有效避免漏掉真正的門限值。3.3 固定滯后階數(shù)與門限延遲用 AIC/BIC 在模型空間里選擇門限值定下來之后還需要確定兩件事每個體制內(nèi)的 AR 滯后階數(shù) p以及門限延遲 d如果用序列自身做門限變量就是 z_t y_{t-d} 里的 d。這兩組參數(shù)的最佳組合通??啃畔?zhǔn)則搜索確定。常見的搜索范圍是 p ∈ 1~5d ∈ 1~3組合數(shù)量不大暴力枚舉就行。import warnings warnings.filterwarnings(ignore) def select_tar_order(series, p_range, d_range, threshold_candidates, min_samples20): 網(wǎng)格搜索使 AIC 最小的 (p, d) 組合 results [] for p, d in itertools.product(p_range, d_range): # 門限變量使用 y_{t-d}序列做回歸時使用 1..p 期滯后 try: res fit_tar_rss(series, delayd, threshold_candidatesthreshold_candidates, ar_orderp) if res[threshold] is None: continue model_low, model_high, _, _ res[models] aic_total model_low.aic model_high.aic results.append((aic_total, p, d, res[threshold], res[rss])) except Exception: continue results.sort(keylambda x: x[0]) best_aic, best_p, best_d, best_th, best_rss results[0] print(f最優(yōu) AIC{best_aic:.2f}, p{best_p}, d{best_d}, 門限{best_th:.4f}) return results result_sorted select_tar_order( series_clean, p_rangerange(1, 6), d_rangerange(1, 4), threshold_candidatesseries_clean.shift(1).dropna().quantile(np.arange(0.2, 0.8, 0.1)).values )這里我用AIC而不是BIC—— TAR 模型本身分段估計參數(shù)已經(jīng)比線性模型多了一倍再用 BIC 這種偏保守的準(zhǔn)則容易選出過簡模型AIC 在樣本量中等的情況下更兼顧擬合與復(fù)雜度。不過如果你的樣本量很大超過 500 個觀測換成 BIC 問題也不大最終效果差異很小。還有個細(xì)節(jié)AutoReg的ssr屬性和aic屬性直接獲取殘差平方和與 Akaike 信息準(zhǔn)則值省去手寫公式的麻煩。跑完網(wǎng)格后注意檢查最優(yōu)組合和第二優(yōu)組合的 AIC 差距如果差距很小小于 2說明模型選擇不夠穩(wěn)健需要檢查數(shù)據(jù)里有沒有異常值擾動。4. jasa_03m 月度數(shù)據(jù)的實戰(zhàn)建模從檢驗到預(yù)測的完整流程4.1 為什么月度數(shù)據(jù)特別適合門限模型我拿 jasa_03m 這批月度數(shù)據(jù)舉例因為它非常有代表性。月度觀測天然具有三種結(jié)構(gòu)第一是季節(jié)性周期比如 12 個月的景氣循環(huán)第二是趨勢變化比如經(jīng)濟擴張期和收縮期第三是變點事件比如疫情沖擊、政策調(diào)整、供需關(guān)系的結(jié)構(gòu)性轉(zhuǎn)變。這三者疊加在一個線性模型里往往互相掩蓋——季節(jié)項把突變點抹平趨勢項把短期的門限效應(yīng)吞掉。門限模型處理月度數(shù)據(jù)的優(yōu)勢在于它能自動根據(jù)門限變量可以是滯后值也可以是外部月度指標(biāo)把數(shù)據(jù)分成“高體制”和“低體制”比如物價高企的月份和物價低迷的月份兩個體制內(nèi)部分別擬合不同的自回歸結(jié)構(gòu)捕捉到的依賴關(guān)系就完全不一樣。你可能已經(jīng)注意到這有點像把“時間分段”換成“狀態(tài)分段”——不是按時間軸切一刀而是按變量取值切一刀這樣即便相同月份處于不同年份只要狀態(tài)相似就會被分到同一個體制里估計參數(shù)信息的利用率高得多。用 jasa_03m 這類數(shù)據(jù)做實證時我通常先把序列拆成訓(xùn)練集和測試集比如前 80% 訓(xùn)練后 20% 驗證在訓(xùn)練集上完成全部參數(shù)選擇測試集只用于最終的滾動預(yù)測評估。4.2 帶外生門限變量的 TAR把外部驅(qū)動因素引入模型做月度數(shù)據(jù)時門限變量只選序列自身的滯后值有時不夠。比如你預(yù)測一個城市的月度用電量門限變量用“上月氣溫是否超過某個閾值”比用“上個月用電量是否超過某個值”更貼近物理實際。這就回到 2.2 節(jié)講的 TAR 和 SETAR 的分野。帶外生門限變量的實現(xiàn)方式與純 SETAR 幾乎一樣只是門限變量的數(shù)據(jù)來源變了。def fit_tar_exog(series, exog_threshold, threshold_candidates, ar_order): 外生門限變量的 TAR 擬合 exog_threshold: 與 series 等長的外生門限變量如氣溫、利率 best {rss: np.inf, threshold: None, models: None} for th in threshold_candidates: mask_low exog_threshold th mask_high exog_threshold th s_low, s_high series.loc[mask_low].dropna(), series.loc[mask_high].dropna() if len(s_low) 20 or len(s_high) 20: continue model_low AutoReg(s_low, lagsar_order).fit() model_high AutoReg(s_high, lagsar_order).fit() rss_total model_low.ssr model_high.ssr if rss_total best[rss]: best[rss] rss_total best[threshold] th best[models] (model_low, model_high, mask_low, mask_high) return best和 3.2 的fit_tar_rss一比唯一的區(qū)別是第四行的exog_threshold替代了原有的series.shift(delay)其他一切照舊??雌饋砀膭硬淮蟮@個選擇的背后是建模思路的變化當(dāng)門限變量來自外部時你得先確認(rèn)它的“時間對齊”——是同步值當(dāng)月氣溫還是滯后值上月氣溫這取決于業(yè)務(wù)上哪個變量“驅(qū)動”了序列切換機制。做外生門限變量時最容易翻車的就是對齊問題門限變量和序列的索引對應(yīng)錯位一位整個模型就串味了。我在項目里會先用print(pd.DataFrame({series: series, gate: exog_threshold}).dropna().head())檢查對齊再做擬合這一步五分鐘能省三小時的排查功夫。4.3 模型預(yù)測與雙體制的樣本外對比模型估計完最終要回答的問題是它比單一線性 AR 模型強在哪兒驗證方法不復(fù)雜——用訓(xùn)練好的門限模型對測試集逐期做預(yù)測并和線性 AR 模型做誤差對比。門限模型的預(yù)測邏輯是每一步預(yù)測前先看門限變量落在哪個體制再用對應(yīng)體制的 AR 模型來預(yù)測。預(yù)測值可以直接用AutoReg的predict方法但門限模型的預(yù)測要小心滯后階數(shù)的對齊推薦把預(yù)測函數(shù)封裝好。def tar_forecast(series, exog_threshold, th, model_low, model_high, ar_order, horizon): 門限模型滾動預(yù)測, h 步預(yù)測 forecasts [] history list(series.values) gate_history list(exog_threshold.values) for _ in range(horizon): gate_val gate_history[-1] # 當(dāng)前狀態(tài)決定用哪個體制的模型預(yù)測 active_model model_low if gate_val th else model_high # 取最近 ar_order 個觀測作為起點 recent history[-ar_order:] # 手動構(gòu)造滯后特征 X np.array([recent[::-1]]).reshape(1, -1) yhat active_model.predict(X)[0] if hasattr(active_model, predict) else np.nan # 手動用參數(shù)算滯后系數(shù) 常數(shù)項 if np.isnan(yhat): yhat active_model.params.iloc[0] # 常數(shù) for lag_i in range(1, ar_order1): yhat active_model.params.iloc[lag_i] * recent[-lag_i] forecasts.append(yhat) # 更新歷史序列 history.append(yhat) gate_history.append(gate_val) return forecasts這段代碼用的active_model.predict(X)是理論化寫法實際用statsmodels的AutoReg對象做單步預(yù)測時需要先擬合一個AutoReg的預(yù)測包裝器或者直接讀取參數(shù)做滾動推算。我更推薦后一種方式把active_model.params按常數(shù)和滯后系數(shù)拆開手動加權(quán)計算下一期預(yù)測值邏輯透明、不容易因為時序索引錯位而出錯。horizon是預(yù)測步長月度數(shù)據(jù)一般做 3~12 步預(yù)測。預(yù)測結(jié)束后對比mean_absolute_error(forecasts, actual)和線性 AR 的同一指標(biāo)如果門限模型沒有顯著優(yōu)勢就回頭檢查門限值是不是不夠顯著或者門限體制的樣本量差異太大——這些在下一章詳細(xì)說。5. 門限自回歸避坑指南五個最容易翻車的地雷5.1 現(xiàn)象門限效應(yīng)檢驗的 p 值不穩(wěn)定換個樣本區(qū)間結(jié)果就變原因門限效應(yīng)檢驗依賴自舉抽樣的隨機性不同的隨機種子會產(chǎn)生不同的模擬分布p 值落在 0.04 到 0.06 之間浮動時結(jié)論就會搖擺。另一個原因是樣本量不夠自舉分布的形狀不夠穩(wěn)定。解決跑檢驗時固定隨機種子np.random.seed(42)并且至少跑三次不同的種子觀察 p 值是否穩(wěn)定地在 0.05 以下。如果 p 值在臨界點附近晃動別急著建模先用數(shù)據(jù)可視化確認(rèn)門限變量的分布是否存在雙峰或明顯的跳躍結(jié)構(gòu)。我在月度數(shù)據(jù)上踩過這個坑第一次跑 p 值是 0.03換了個子樣本變成 0.08最后發(fā)現(xiàn)是樣本里有一個極端的異常值把門限點拉偏了剔除后 p 值穩(wěn)定在 0.01。5.2 現(xiàn)象模型在訓(xùn)練集上擬合得很好測試集上一塌糊涂原因這是門限模型最典型的過擬合陷阱——搜索門限值的時候用的是“殘差平方和最小”這個目標(biāo)但是候選門限點本身是從觀測值里挑的你其實是在用測試信息反選門限。本質(zhì)上你是在拿數(shù)據(jù)來擬合門限的位置而不是先驗地指定它。這會導(dǎo)致門限點選擇過度優(yōu)化訓(xùn)練集內(nèi)部再分段做回歸自由度大幅上升泛化能力自然大幅下降。解決這一點我算交過學(xué)費的。后來固定做法是——門限值只從訓(xùn)練集的前 80% 里選后 20% 完全不參與門限值搜索選完門限后再用全訓(xùn)練集估參數(shù)。這樣至少保證門限點的選擇沒有“偷看”最后一段信息。另外給每個體制設(shè)定最小觀測數(shù)我通常設(shè)為總樣本的 15% 以上防止門限切在邊緣導(dǎo)致單體制樣本太少。5.3 現(xiàn)象兩個體制的滯后階數(shù)都用 AIC 定但模型整體比線性 AR 還差原因每個體制單獨用 AIC 選滯后階數(shù)可能在低體制選出 p3高體制選出 p5兩套模型的復(fù)雜度加在一起遠(yuǎn)超數(shù)據(jù)能支撐的信息量造成過擬合。更隱蔽的問題是門限模型每個體制內(nèi)的觀測數(shù)少相同階數(shù)的 AR 模型在小樣本下估計方差更大。解決我一般會先固定兩個體制使用相同的滯后階數(shù)再用聯(lián)合 AIC兩個體制 AIC 之和去選。如果不同體制的滯后階數(shù)差異在業(yè)務(wù)上站不住腳就保持同構(gòu)。還有一種更實用的策略殘差診斷用 Ljung-Box 檢驗如果殘差已經(jīng)近似白噪聲就不要再加滯后階數(shù)了——AIC 在這個場景下容易鉤住多余的主項。5.4 現(xiàn)象門限值估計出來正好落在數(shù)據(jù)的中位數(shù)附近但切換邏輯說不通原因這不一定是你選錯了但也可能是門限變量選錯了。比如用 y_{t-1} 做門限門限值切出來的兩段剛好把數(shù)據(jù)按高低排序切開但這兩段沒有業(yè)務(wù)上的“機制差異”只是數(shù)學(xué)上的偶然。這種模型的預(yù)測能力往往很脆弱換一個數(shù)據(jù)集門限值就漂移。解決回到業(yè)務(wù)層面問一個問題這個門限變量的哪個數(shù)值在你的領(lǐng)域里意味著“機制切換”比如對銷量數(shù)據(jù)門限可能是“上期銷量高于歷史均值 1.2 倍”——這是產(chǎn)能瓶頸點對風(fēng)速數(shù)據(jù)門限可能是“切變層高度超過某個海拔”——這是風(fēng)切變機制切換點。如果門限值給不出業(yè)務(wù)解釋就換門限變量別硬撐著用數(shù)學(xué)結(jié)果。5.5 現(xiàn)象帶外生門限變量時兩種體制的參數(shù)看起來沒有顯著差異原因這說明門限變量選的沒有區(qū)分度——它在兩個區(qū)間內(nèi)對應(yīng)的序列行為本質(zhì)上一樣分段只是給模型額外塞了一堆參數(shù)。另一個常見原因是外生門限變量和序列同期相關(guān)導(dǎo)致門限效應(yīng)被序列自身的慣性吸收門限切換本身沒有增量信息。解決在擬合門限變量之前先分別統(tǒng)計兩個體制內(nèi)序列的均值、方差和自相關(guān)系數(shù)。如果三個統(tǒng)計量沒有顯著差異就不要做門限模型??梢杂煤唵蔚?t 檢驗比較兩個體制的均值差異是否顯著。我自己有個鐵律門限模型如果只是讓 AIC 降低了不到 5%這個模型就沒有實用價值不如老老實實用線性 ARIMA。6. 殘差自舉法計算預(yù)測區(qū)間門限模型的不確定性量化最后一節(jié)是門限模型里最實用但很多人沒做到位的一環(huán)預(yù)測區(qū)間。區(qū)間的價值在于幫你區(qū)分“正常的波動”和“異常的偏移”。在月度數(shù)據(jù)上一個未來 3 個月的預(yù)測值如果沒有區(qū)間業(yè)務(wù)方根本不知道數(shù)字是大概率事件還是賭運氣。殘差自舉法的操作分四步第一步從門限模型的殘差序列里隨機抽樣第二步把抽樣殘差疊加到預(yù)測值上生成一條模擬的未來路徑第三步重復(fù) 1000 次收集所有模擬路徑在每個未來時間點的值第四步取 2.5% 和 97.5% 分位數(shù)作為 95% 置信區(qū)間。這個方法的優(yōu)勢是不需要假設(shè)殘差服從正態(tài)分布——月度數(shù)據(jù)經(jīng)常有偏態(tài)用正態(tài)假設(shè)會算出過窄的區(qū)間。def tar_forecast_interval(residuals, series, exog_threshold, th, model_low, model_high, ar_order, horizon, n_sim1000): 殘差自舉法生成預(yù)測區(qū)間 sim_paths [] for _ in range(n_sim): # 有放回抽樣殘差 boot_resid np.random.choice(residuals, sizehorizon, replaceTrue) # 用不帶噪聲的模型預(yù)測初始路徑 base_forecast tar_forecast(series, exog_threshold, th, model_low, model_high, ar_order, horizon) # 疊加殘差生成模擬路徑 sim_path np.array(base_forecast) boot_resid sim_paths.append(sim_path) sim_paths np.array(sim_paths) lower np.percentile(sim_paths, 2.5, axis0) upper np.percentile(sim_paths, 97.5, axis0) return np.array(base_forecast), lower, upper代碼里最關(guān)鍵的一行是np.random.choice(residuals, sizehorizon, replaceTrue)它決定了重抽樣的隨機性來源。如果殘差序列存在自相關(guān)直接重抽樣會破壞時序結(jié)構(gòu)此時應(yīng)該先對殘差擬合一個 AR(1) 模型再抽取白噪聲部分做自舉。我在月度數(shù)據(jù)上的經(jīng)驗是殘差如果有明顯的序列相關(guān)Ljung-Box 檢驗 p 值小于 0.05先對殘差做 ARIMA 過濾再對過濾后的白噪聲自舉區(qū)間會準(zhǔn)很多。預(yù)測區(qū)間的寬度是模型不確定性的直觀標(biāo)尺如果區(qū)間寬度超過預(yù)測值本身的 50%說明模型在這個數(shù)據(jù)集上的信息量有限后面要考慮引入外部變量或更高的采樣頻率。個人習(xí)慣是任何門限模型的預(yù)測報告里一律把單點預(yù)測和區(qū)間一起交付。區(qū)間不是錦上添花是判斷模型可靠性的核心依據(jù)——如果一位工程師只給我單點預(yù)測卻不給區(qū)間我基本不會采用。希望這些內(nèi)容能幫你在做門限自回歸的路上少走些彎路把這套非線性工具真正用起來。本文還有配套的精品資源點擊獲取