此文內(nèi)容為《動態(tài)隨機(jī)一般均衡(DSGE)模型》的筆記浪听，李向陽老師著，清華大學(xué)出版社出版。
我只將個人會用到的知識作了筆記宫仗，并對教材較難理解的部分做了進(jìn)一步闡述。為了更易于理解旁仿，我還對教材上的一些部分（包括代碼和正文）做了修改藕夫。
僅供學(xué)習(xí)參考，請勿轉(zhuǎn)載枯冈，侵刪毅贮！
上一篇有關(guān)Dyanre的文章為大家介紹了有關(guān)確定性模擬的技術(shù)。在本次推文尘奏，我將介紹隨機(jī)模擬技術(shù)滩褥。另外，我還會介紹有關(guān)參數(shù)估計的問題炫加。
3.9 隨機(jī)模擬分析：stock_simul
隨機(jī)模擬是Dyanre最常用的計算命令之一瑰煎。隨機(jī)模擬首先使用擾動算法完成模型的求解任務(wù)铺然，然后再求解的基礎(chǔ)上計算內(nèi)生變量的脈沖響應(yīng)、各階矩丢间，最后進(jìn)行隨機(jī)模擬探熔，——即根據(jù)隨機(jī)抽取的外生沖擊的樣本，模擬出內(nèi)生變量的增長路徑烘挫。
3.9.1 隨機(jī)模擬命令簡介
在Dynare中诀艰，隨機(jī)模擬的選項多達(dá)44個，常用的有：
drop = xx饮六，計算內(nèi)生變量的各階矩時其垄，需要首先指定丟棄的期數(shù)，默認(rèn)為100hp_filter = xx卤橄，計算內(nèi)生變量的各階矩時绿满，首先使用HP濾波，xx代表值
irf=xx窟扑，計算xx期的脈沖響應(yīng)喇颁，默認(rèn)為40irf_shocks=(shock1, shock2,...)，針對給定的外生沖擊計算脈沖響應(yīng)嚎货，默認(rèn)計算全部irf_plot_threshold=xx橘霎，指定畫脈沖響應(yīng)圖形的門檻大小，默認(rèn)nocorr殖属，指定不要在Matlab窗口顯示相關(guān)系數(shù)矩陣nofunctions姐叁，指定不要在Matlab窗口顯示政策函數(shù)和轉(zhuǎn)換方程nomoments，指定不要在Matlab命令窗口顯示各階矩nograph洗显，指定不要顯示和保存圖形nodisplay外潜，指定不要顯示圖形，但保存到硬盤noprint挠唆，指定不要顯示任何計算信息处窥，對循環(huán)很有幫助order=xx，指定Taylor近似對階數(shù)玄组，可為1``2``3滔驾，默認(rèn)為2periods=xx，指定隨機(jī)模擬的期數(shù)巧勤，同時也是決定計算理論矩還是模擬矩的關(guān)鍵所在pruning嵌灰，指定剪枝算法弄匕，即在迭代過程中丟棄高階項
3.9.2 隨機(jī)模擬示例
以第三章第一節(jié)的RBC模型為例颅悉，這里示例如何使用模擬結(jié)果中的水平變量進(jìn)行做圖分析：
/*
* 《動態(tài)隨機(jī)一般均衡(DSGE)模型》教材
* 隨機(jī)模擬示例
* 2020年3月30日
* @愛吃漢堡薯條
* @侵刪
*/
%----- 源代碼9 -----%
var c k lab z I; % 內(nèi)生變量聲明
varexo e; % 外生變量聲明
parameters beta theta delta alpha tau rho s; % 參數(shù)聲明
% 參數(shù)賦值
beta = 0.987;
theta = 0.357;
delta = 0.012;
alpha = 0.4;
tau = 2;
rho = 0.95;
s = 0.01;
%----- 源代碼10 -----%
model;
% 定義模型的“局部變量”，以提高編程的便利性
# mar_c = (c^theta * (1-lab)^(1-theta) ) ^ (1-tau); % 局部變量的聲明使用“#”命令
# mar_c1 = (c(+1)^theta * (1-lab(+1))^(1-theta) ) ^ (1-tau);
% (1) 歐拉方程
[name = 'Euler equation'] % 這是系統(tǒng)標(biāo)簽(tag)迁匠，可以對每個均衡條件進(jìn)行標(biāo)識
mar_c/c = beta * (mar_c1/c(+1)) * (1 + alpha*exp(z)*k^(alpha-1)*lab(+1)^(1-alpha) - delta);
% (2) 勞動供給方程
[name = 'Labor supply']
((1-theta)/theta) * (c/(1-lab)) = (1-alpha)*exp(z)*k(-1)^alpha*lab^(-alpha);
% (3) 資源約束方程
[name = 'Resource constraint']
c + k - (1-delta)*k(-1) = exp(z)*k(-1)^alpha*lab^(1-alpha);
% (4) 資本積累方程
[name = 'Capital accumulation']
I = k - (1-delta)*k(-1);
% (5) 外生沖擊AR(1)過程
[name = 'Technology shock']
z = rho*z(-1) + s*e;
end;
%----- 源代碼11 -----%
initval; % 初值模塊剩瓶，設(shè)定內(nèi)生變量和外生變量的初始值
k = 0.5;
c = 0.5;
lab = 0.5;
z = 0;
I = 0.5;
e = 0;
end;
shocks; % 外生沖擊設(shè)定驹溃，直接設(shè)定外生沖擊的方差
var e = s^2;
end;
steady; % 計算穩(wěn)態(tài)
resid; % 給定內(nèi)生變量的取值（通常為穩(wěn)態(tài)）唆姐，計算均衡條件對應(yīng)的靜態(tài)方程的殘差
%----- 源代碼12 -----%
% 如果不指定模擬的次數(shù)“periods”僵娃，Dynare就不會進(jìn)行模擬
stoch_simul(periods=1000, irf=400, order=1);
% 將模擬結(jié)果保存為“simudata.mat”文件
dynasave('simudata.mat');
則有下圖：
在這里，隨機(jī)模擬命令指定了3個選項：
模擬1000期
計算40期脈沖響應(yīng)圖
進(jìn)行二階泰勒近似求解
前面已經(jīng)提到鸭廷，內(nèi)生變量的模擬路徑結(jié)果儲存于oo_.endo_simul結(jié)構(gòu)數(shù)組中枝缔。其行的順序?qū)?yīng)著以聲明順序排序的內(nèi)生變量布疙，其列則為模擬期數(shù)。
其中愿卸，技術(shù)沖擊的均值為：
>> mean(oo_.endo_simul(4,:))
ans = -7.2640e-05
3.10 參數(shù)估計簡介
模型中結(jié)構(gòu)參數(shù)的賦值方法有兩種灵临，分別是：
- 校準(zhǔn)（Calibration），直接賦值
- 估計（Estimation）趴荸，使用實際統(tǒng)計或觀測數(shù)據(jù)進(jìn)行估計
簡單來說儒溉，Dynare對參數(shù)估計就是使用某些方法，比如貝葉斯估計发钝、極大似然估計顿涣，從數(shù)據(jù)中發(fā)現(xiàn)真實取值
3.10.1 極大似然估計與貝葉斯估計的基本邏輯
本文僅簡單闡述極大似然估計和貝葉斯估計的基本原理和邏輯，對背后的復(fù)雜數(shù)理知識酝豪，請參考教材和博客涛碑。
3.10.1.1 極大似然估計
極大似然估計（Maximum Likelihood Estimation, MLE）就是選擇參數(shù) 來最大化如下的似然函數(shù)（Likelihood），其原理比較簡單：
其中寓调，條件概率 是使用
濾波 進(jìn)行計算的锌唾。
總而言之，它是選擇 使已經(jīng)發(fā)生時間的概率最大夺英。
3.10.1.2 貝葉斯估計
貝葉斯估計（Bayesian Estimation）的基本原理非常直觀晌涕，就是在參數(shù)先驗分布（Prior）的基礎(chǔ)上，結(jié)合數(shù)據(jù)的信息（）痛悯，找到參數(shù)的后驗分布（Posterior）余黎。也就是說，貝葉斯估計是一種簡單的
映射载萌，依據(jù)數(shù)據(jù)信息惧财，將先驗分布映射為后驗分布。
由于先驗分布是依據(jù)經(jīng)驗而假設(shè)的分布扭仁，可能存在一定的偏誤垮衷，而數(shù)據(jù)的作用其實是將數(shù)據(jù)中有用的信息融入，以糾正偏誤并形成后驗分布乖坠。
從邏輯上說搀突，貝葉斯估計的背后就是貝葉斯法則，即參數(shù)的后驗分布由條件概率公式：
推導(dǎo)出：
從數(shù)據(jù)使用和先驗假設(shè)的角度看熊泵，貝葉斯估計位于極大似然估計和參數(shù)校準(zhǔn)之間：
- 極大似然估計：直接使用數(shù)據(jù)中的信息仰迁，沒有先驗
- 貝葉斯估計：即使用數(shù)據(jù)中的信息甸昏，也具有先驗信息
- 校準(zhǔn)：僅使用先驗信息，沒有使用數(shù)據(jù)的信息
貝葉斯估計的先驗分布是一般意義上的分布徐许，即標(biāo)準(zhǔn)差不為0施蜜，不是一個點，而是一個區(qū)域雌隅。嚴(yán)格地說翻默，是整個參數(shù)空間的一個子區(qū)域，作為候選區(qū)域恰起。先驗分布可以理解為施加于似然函數(shù)的“權(quán)重”或者“偏好”冰蘑，使得參數(shù)的選擇范圍集中于均值附近。
換句話說村缸，先驗分布使得“巨大的參數(shù)選擇空間”變小祠肥，使似然函數(shù)集中于“更有意義”的子區(qū)域去尋找合適的參數(shù)估計值。Dynare的貝葉斯分布可分為兩步：
- 使用數(shù)值求解算法梯皿，計算似然函數(shù)取最大值的
眾數(shù)（Model） - 從
眾數(shù)開始仇箱，使用MCMC方法模擬后驗分布，然后計算各種感興趣的矩（Moments）
3.10.2 馬爾可夫鏈——蒙特卡洛（MCMC）算法
由于參數(shù)的后驗分布（密度函數(shù)）在大多數(shù)DSGE模型中（甚至某些RBC模型中）無法使用解析式的形式來表達(dá)东羹，因此不能通過貝葉斯法則計算出來剂桥，只能通過隨機(jī)抽樣的方法，即蒙特卡洛算法（Monte Carlo）來近似逼近后驗分布属提。
馬爾可夫鏈-蒙特卡洛方法（Markov Chain-Monte Carlo, MCMC）是一類方法的統(tǒng)稱权逗，是貝葉斯估計的邏輯方法和框架，非具體估計的算法冤议≌遛保可以參考 知乎 網(wǎng)友的文章。實際上我們不太需要了解背后的數(shù)學(xué)原理恕酸。
關(guān)于本節(jié)堪滨，我們最重要的是基本了解Dynare貝葉斯估計的核心邏輯，如圖3.12所示：
3.10.3 Dynare參數(shù)估計和一個例子
我們著重于介紹如何用Dynare內(nèi)置的估計命令蕊温，對一個或多個參數(shù)進(jìn)行編程估計袱箱。
3.10.3.1 Dynare參數(shù)估計的基本邏輯
Dynare參數(shù)估計是使用以下模型的一階線性化系統(tǒng)進(jìn)行估計：
其中， 是待估計參數(shù)向量义矛；
是外生沖擊发笔；
是內(nèi)生變量組成的向量。
Dynare求解參數(shù)估計的基本邏輯是：
- 計算內(nèi)生變量
的穩(wěn)態(tài)值
- 線性化系統(tǒng)均衡條件
- 求解線性化模型
- 使用Kalman濾波算法計算對數(shù)似然函數(shù)（Log-Likelihood）
- 找到極大似然值分布的眾數(shù)（Model）
- 使用Metropolis-Hastings算法模擬后驗分布
3.10.3.2 Dynare參數(shù)估計的語法
Dyanre參數(shù)估計的基本語法分成三部分：
- 聲明參數(shù)估計模塊
- 待估參數(shù)本身的聲明
- 使用
estimation命令
Dynare中聲明估計參數(shù)模塊的內(nèi)置命令是以estimated_params;開頭凉翻，以end;結(jié)束了讨，中間輸入待估計參數(shù)聲明。待估計參數(shù)有兩種常見的類型：
- 結(jié)構(gòu)參數(shù)
- 外生參數(shù)的標(biāo)準(zhǔn)差參數(shù)
具體實例如下：
estimated_params; % 結(jié)構(gòu)參數(shù)的估計
gamma, normal_pdf, 1,0.05;
alpha, uniform_pdf, , , 0, 1; % uniform分布需要留逗號，見下文
end;
estimated_params; % 外生沖擊標(biāo)準(zhǔn)差參數(shù)的估計
stderr e, inv_gamma_pdf, 0.01, inf;
end;
結(jié)構(gòu)參數(shù)的聲明是：
- 以待估計參數(shù)的名稱開頭
- 然后是該參數(shù)服從的先驗分布的名稱量蕊，該名稱由Dynare預(yù)先定義，如表3.11所示：
| 先驗分布名稱 | 符號形式 | 分布取值范圍 |
|---|---|---|
beta_pdf |
||
gamma_pdf |
||
inv_gamma_pdf |
||
normal_pdf |
||
uniform_pdf |
其中艇挨， 和
分別表示均值和標(biāo)準(zhǔn)差残炮；參數(shù)
、
僅用于
Beta分布缩滨、 Gamma分布和Uniform分布势就，可以為無窮大。此外脉漏，Uniform分布因為不需要聲明均值和標(biāo)準(zhǔn)差苞冯，因此需要在參數(shù)聲明中留出兩個逗號隔開的空格，然后聲明 侧巨、
參數(shù)舅锄。
標(biāo)準(zhǔn)差參數(shù)的聲明：
- 需要使用關(guān)鍵字
stderr后跟外生沖擊的名稱，而不是外生沖擊對應(yīng)的標(biāo)準(zhǔn)差的名稱 - 然后是指定其服從的先驗分布司忱，聲明語法與前面的相同
3.10.3.2 注意事項
對于先驗分布的選擇皇忿，確實是一個困難的問題。選擇先驗分布具有主觀性坦仍，這成為貝葉斯分布受質(zhì)疑和批評的主要原因鳍烁。一般參考：
- 對于標(biāo)準(zhǔn)差
的估計，通常選擇
Inverse Gamma 分布 - 對于
過程的滯后參數(shù)
通常選擇
Beta分布 - 其他參數(shù)繁扎，大多使用
Normal分布
在估計之前幔荒，還要注意以下幾個問題：
數(shù)據(jù)問題：估計參數(shù)的前提是模型中某些內(nèi)生變量是可以觀測到或者由觀測值經(jīng)過處理得到，且觀測變量的個數(shù)不能超過模型中外生沖擊的個數(shù)梳玫，否則無法計算模型的似然函數(shù)爹梁。
觀測數(shù)據(jù)與內(nèi)生變量的匹配：在估計時，要求模型文件的變量形式和觀測數(shù)據(jù)相匹配提澎。
-
觀測變量的聲明問題：內(nèi)生變量如果聲明為可觀測變量卫键，則需要使用以下命令
% 語法 varobs [VARIABLE_NAME]; % 例如 varobs pi tau a; % 聲明三個內(nèi)生變量為觀測變量
3.10.3.3 estimation命令介紹
Dynare會在estimation命令中調(diào)取外部數(shù)據(jù)（使用datafile選項指定外部數(shù)據(jù)文件），并加載到內(nèi)存中去虱朵。觀測變量聲明命令varobs會告訴Dynare去尋找變量名為pi莉炉、tau、a的三個序列碴犬，如果找不到絮宁，會提示錯誤。
Estimation`命令及其選項問題：這里僅介紹幾個最重要的關(guān)鍵選項
-
datafile=xxx：可選擇Matlab的.m文件服协、.mat矩陣或者.csv文件等 -
conf_sig=xxx：定義置信區(qū)間的寬度（如xxx=.95定義了“95%的置信區(qū)間”） -
frst_obs=xxx：定義由load_Gali_est_data.m加載的數(shù)據(jù)绍昂，第一個開始用于對估計的數(shù)據(jù)位置，xxx必須是一個正整數(shù)，缺省為1 -
nobs=yy：Dynare估計中使用的觀測變量序列長度（t=xxx到t=yy+xxx-1窘游，同樣的yy必須是正整數(shù)唠椭，缺省為全部觀測變量 -
mode_check：畫出后驗分布眾數(shù)附近的后驗分布密度圖 -
mode_compute=xxx：選擇優(yōu)化算法，取值0~10忍饰，缺省使用4（該命令在不斷更新中） -
forecast=n：在觀測變量后一期之后預(yù)測期值
-
smmother：計算內(nèi)生變量和外生沖擊的后驗分布的相關(guān)矩 -
mh_replic=xxx：定義計算后驗分布時MCMC重復(fù)的次數(shù)贪嫂，缺省為2000 -
mh_jscale=xxx：定義MCMC中的乘法因子。理想情況下應(yīng)該選擇使“接受率”落在25%~33%之間艾蓝。該參數(shù)需要多次調(diào)試后獲得比較合適的值力崇，缺省為0.2 -
mh_nblocks=xxx：定義Metropolis-Hastings算法平行進(jìn)行的馬爾可夫鏈的個數(shù)，缺省為2 -
filtered_vars：計算濾波變量赢织，結(jié)果存儲在oo_.FilteredVarialbes中 -
moments_varendo：計算內(nèi)生變量的后驗分布的理論矩亮靴，結(jié)果存儲在oo_.PosteriorTheoreticalMoments中
3.10.3.4 極大似然估計的例子
考慮一個新凱恩斯模型，不含投資但含粘性價格和產(chǎn)出缺口于置。
（1）模型
其中：
-
為消費跨期替代彈性的倒數(shù)
-
為勞動負(fù)效用的外生沖擊茧吊，服從
-
為勞動供給的
彈性的倒數(shù)
首先是新凱恩斯IS曲線（NKIS）和新凱恩斯菲利普斯曲線（NKPC）：
其次是自然水平產(chǎn)出和自然利率方程：
假設(shè)貨幣沖擊政策為Taylor規(guī)則：
從而，模型由 和六個均衡條件：
- NKIS曲線
- NKPC曲線
- 自然利率方程
- Taylor規(guī)則
- 技術(shù)沖擊的
- 偏好沖擊的
接下來八毯，對上述模型文件進(jìn)行隨機(jī)模擬饱狂，模擬5000個樣本，作為下一步估計的數(shù)據(jù)宪彩。然后使用上述模擬的樣本作為觀測數(shù)據(jù)休讳，對模型中3個結(jié)構(gòu)參數(shù)進(jìn)行估計，如表3.12所示尿孔，以示例如何在Dynare中分別使用極大似然估計和貝葉斯估計估計此三個參數(shù)：
| 待估計參數(shù) | 含義 |
|---|---|
| 技術(shù)沖擊 |
|
| 勞動負(fù)效用外生沖擊 |
|
| 技術(shù)沖擊 |
隨機(jī)模擬的代碼如下（請首先運行）：
// 生成用于參數(shù)估計的隨你模擬序列
//Written by Xiangyang Li@SCC
//declaration of the endogenous variables, log-devation form
var a // technology
pi //CPI inflation
i //nominal rate
istar //nature real rate
tau //labor disutility shock
x //output gap = output - natural output
dy //growth rate of output
;
varexo eps_a eps_tau eps_i;
parameters
sigma phi beta kappa phi_x phi_pi rho lambda rhoi theta;
// Parameter Values
sigma = 1; //inverse of inter-temporal elasticity of substitution of consumption
beta = 0.97; //discount factor of households
phi_x = .15; //coefficient of output gap in monetary policy rule
phi_pi = 1.5; //coefficient of inflation in monetary policy rule
rhoi = 0.8; //persistence of nominal interest rate;
rho = 0.8; //persistence of technology
lambda = 0.5; //persistence of labor disutility shock;
phi = 1; //inverse of Frisch elasticity of labor supply
theta = 0.75; //Calvo stickiness parameter
kappa = ((1-theta)*(1-beta*theta))/theta; //simplifying parameter
// The model in log-linear
model(linear);
//(1) Philiphs Curve - Calvo Pricing Equation
beta*pi(+1) + kappa*x = pi;
//(2) New Keynesian IS curve
sigma*(i - pi(+1)-istar) = x(+1) - x;
//(3) natural real rate definition
istar = sigma*(1+phi)*(rho-1)/(sigma+phi)*a + sigma*(1-lambda)/(sigma+phi)*tau;
//(4)Taylor Rule
i= rhoi*i(-1) + (1-rhoi)*(phi_pi*pi + phi_x*x) +eps_i;
//(5)Technology shock
a = rho*a(-1) + eps_a;
//(6) Labor disutility shock
tau = lambda*tau(-1) + eps_tau;
//(7) Growth rate of output - by output gap
dy = x - x(-1) + (1+phi)/(sigma+phi)*(a-a(-1)) - (tau - tau(-1))/(sigma+phi); // construction of output growth
end;
shocks;
var eps_a;
stderr .02;
var eps_tau;
stderr .02;
var eps_i;
stderr 0.02;
end;
set_dynare_seed=1;
stoch_simul(periods=5000, irf=7, nograph);
gcsim = oo_.endo_simul;
save data gcsim;
（2）最大似然估計的實踐
在下面的代碼中俊柔，使用了estimated_param;設(shè)定了初始值，使用estimated_param_bounds;設(shè)定了上下邊界活合，此外使用varobs定義了3個觀測變量（請在運行3.10.3.4的（1）后運行）：
// Maximum Likelihood estimation of the Gali-Christiano Model
//we are going to estimate two persistence parameters, rho and lambda
// and one standard devation parameters
//Written by Xiangyang Li@SCC@2016.12
//Modified by Xiangyang Li@2017.10.12
//declaration of the endogenous variables, log-devation form
var a // technology
pi //CPI inflation
i //nominal rate
istar //nature real rate
tau //labor disutility shock
x //output gap = output - natural output
dy //growth rate of output
;
varexo eps_a eps_tau eps_i;
parameters sigma phi beta kappa phi_x phi_pi rhoi theta;
parameters rho lambda; //to be estimated
// Parameter Values
sigma = 1; //inverse of inter-temporal elasticity of substitution of consumption
beta = 0.99; //discount factor of households
phi_x = .15; //coefficient of output gap in monetary policy rule
phi_pi = 1.5; //coefficient of inflation in monetary policy rule
rhoi = 0.8; //persistence of nominal interest rate;
//rho = 0.8; //persistence of technology, to be estimated
//lambda = 0.5; //persistence of labor disutility shock,to be estimated;
phi = 1; //inverse of Frisch elasticity of labor supply
theta = 0.75; //Calvo stickiness parameter
kappa = ((1-theta)*(1-beta*theta))/theta; //simplifying parameter
// The model in log-linear
model(linear);
//(1) Philiphs Curve - Calvo Pricing Equation
beta*pi(+1) + kappa*x = pi;
//(2) New Keynesian IS curve
sigma*(i - pi(+1)-istar) = x(+1) - x;
//(3) natural real rate definition
istar = sigma*(1+phi)*(rho-1)/(sigma+phi)*a + sigma*(1-lambda)/(sigma+phi)*tau;
//(4)Taylor Rule
i= rhoi*i(-1) + (1-rhoi)*(phi_pi*pi + phi_x*x) +eps_i;
//(5)Technology shock
a = rho*a(-1) + eps_a;
//(6) Labor disutility shock
tau = lambda*tau(-1) + eps_tau;
//(7) Growth rate of output - by output gap
dy = x - x(-1) + (1+phi)/(sigma+phi)*(a-a(-1)) - (tau - tau(-1))/(sigma+phi); // construction of output growth
end;
// 設(shè)置沖擊
shocks;
var eps_tau; stderr 0.01;
var eps_i; stderr 0.01;
end;
// 極大似然估計fa:
// 初始化參數(shù)值雏婶，用于優(yōu)化
estimated_params;
stderr eps_a, 0.01; // 外生沖擊參數(shù)的估計聲明
rho, .80; // 結(jié)構(gòu)參數(shù)的估計聲明
lambda, .50;
end;
// 設(shè)置參數(shù)的邊界
estimated_params_bounds;
stderr eps_a, 0.001, .2;
rho, .001,.95;
lambda, .001,.95;
end;
// 用于估計參數(shù)的觀測源：觀測變量
varobs pi tau a;
// 開始估計
estimation(datafile=load_Gali_est_data,conf_sig =.95,first_obs=101,forecast =46,nobs=4000,mode_check,mode_compute=10) a pi i istar tau x dy;
運行上面的代碼，意味著：
- 使用95%置信區(qū)間
- 從第101個數(shù)據(jù)開始作為觀測值
- 在觀測變量以后預(yù)測46期
- 使用的序列長度是4000期
- 畫出后驗分布眾數(shù)附近的后驗分布密度圖
- 使用第10種優(yōu)化算法（見Manual）
（3）關(guān)于觀測變量的分析
在默認(rèn)情況下白指，極大似然估計的結(jié)果令人滿意留晚，非常接近真實值；但選擇不同的觀測變量告嘲，會影響參數(shù)的點估計错维，如表3.13所示：
直觀上，選擇觀測數(shù)據(jù)時橄唬，要盡量選擇和待估計參數(shù)關(guān)系緊密的數(shù)據(jù)序列赋焕，包括使用盡可能多的數(shù)據(jù)項，以獲得準(zhǔn)確的信息仰楚。
（4）關(guān)于序列長度的分析
同樣隆判，觀測變量序列的長度也會影響估計犬庇。可以看出侨嘀，序列長度縮小100倍時臭挽，標(biāo)準(zhǔn)差約增大十倍。除了參數(shù) 以外咬腕，另外兩個參數(shù)的估計值變化很小欢峰，說明極大似然估計具有一定的穩(wěn)健型：
3.10.3.5 貝葉斯估計的例子
使用貝葉斯分布，要先指定先驗分布郎汪。這里設(shè)定：
- 技術(shù)標(biāo)準(zhǔn)差參
-
- 還設(shè)定了3個參數(shù)的上下界和初始值
如下所示（請在運行3.10.3.4的（1）后運行）：
// Maximum Likelihood estimation of the Gali-Christiano Model
//we are going to estimate two persistence parameters, rho and lambda
// and one standard devation parameters
//Written by Xiangyang Li@SCC@2016.12
//Modified by Xiangyang Li@2017.10.12
//declaration of the endogenous variables, log-devation form
var a // technology
pi //CPI inflation
i //nominal rate
istar //nature real rate
tau //labor disutility shock
x //output gap = output - natural output
dy //growth rate of output
;
varexo eps_a eps_tau eps_i;
parameters sigma phi beta kappa phi_x phi_pi rhoi theta;
parameters rho lambda; //to be estimated
// Parameter Values
sigma = 1; //inverse of inter-temporal elasticity of substitution of consumption
beta = 0.99; //discount factor of households
phi_x = .15; //coefficient of output gap in monetary policy rule
phi_pi = 1.5; //coefficient of inflation in monetary policy rule
rhoi = 0.8; //persistence of nominal interest rate;
//rho = 0.8; //persistence of technology, to be estimated
//lambda = 0.5; //persistence of labor disutility shock,to be estimated;
phi = 1; //inverse of Frisch elasticity of labor supply
theta = 0.75; //Calvo stickiness parameter
kappa = ((1-theta)*(1-beta*theta))/theta; //simplifying parameter
// The model in log-linear
model(linear);
//(1) Philiphs Curve - Calvo Pricing Equation
beta*pi(+1) + kappa*x = pi;
//(2) New Keynesian IS curve
sigma*(i - pi(+1)-istar) = x(+1) - x;
//(3) natural real rate definition
istar = sigma*(1+phi)*(rho-1)/(sigma+phi)*a + sigma*(1-lambda)/(sigma+phi)*tau;
//(4)Taylor Rule
i= rhoi*i(-1) + (1-rhoi)*(phi_pi*pi + phi_x*x) +eps_i;
//(5)Technology shock
a = rho*a(-1) + eps_a;
//(6) Labor disutility shock
tau = lambda*tau(-1) + eps_tau;
//(7) Growth rate of output - by output gap
dy = x - x(-1) + (1+phi)/(sigma+phi)*(a-a(-1)) - (tau - tau(-1))/(sigma+phi); // construction of output growth
end;
shocks;
var eps_tau;
stderr 0.01;
var eps_i;
stderr 0.01;
end;
// 貝葉斯估計
// 設(shè)置先驗分布
estimated_params;
// eps_a的均值為0.02煞赢，標(biāo)準(zhǔn)差為10
stderr eps_a, inv_gamma_pdf, 0.02, 10;
// 設(shè)置另外兩個參數(shù)的先驗
rho, beta_pdf, 0.80, 0.04;
lambda, beta_pdf, 0.80, 0.04;
end;
estimated_params_bounds;
stderr eps_a, 0.001, .2;
rho, .001,.95;
lambda, .001,.95;
end;
// 設(shè)置初值，MLE則不需要這一模塊
estimated_params_init;
stderr eps_a, 0.019;
rho, 0.82;
lambda, 0.51;
end;
// 設(shè)置觀測變量
varobs pi tau a;
estimation(datafile=load_Gali_est_data, // 數(shù)據(jù)文件
conf_sig =.95, // 95%置信區(qū)間
first_obs=101, // 從101位開始使用數(shù)據(jù)
forecast =46, // 在觀測序列后預(yù)測46期
nobs=40, // 使用40期數(shù)據(jù)進(jìn)行預(yù)測
mode_check, // 畫出后驗分布分布密度圖
mode_compute=4, // 使用第4種優(yōu)化算法
mh_replic=1200, // MCMC重復(fù)次數(shù)
mh_jscale=1.4, // MCMC的乘法因子
mh_nblocks=2) // 馬爾可夫鏈個數(shù)
a pi i istar tau x dy;
運行上面的代碼哄孤，可以得到幾個圖照筑，其中下圖：
可以直觀看出貝葉斯估計的先驗和后驗的情況。
貝葉斯估計比極大似然估計更加穩(wěn)绞莩隆（請看教材）凝危，所以一般使用貝葉斯估計。
ng$
