epoh老师，您好！
       我想问一下：fmincon()这个函数是matlab自带的还是jplv7中的函数？
如何调用？

epoh老师，您好！
         您看一下下面的程序运行错在什么地方？我尝试修改，还是不能解决问题。
非常感谢！
          附件中有程序。


>     #######################################################
>     ###       R code for Lec 10 Examples                 ###
>     #######################################################
>     
> ### Example 1 (pepper moth)
>
>
>   moth<-function(p,n.obs){
+     n<-sum(n.obs)
+     nc<-n.obs[1]
+     ni<-n.obs[2]
+     nt<-n.obs[3]
+     ntt<-nt
+   
+     cat(p,"\n")
+     pct<-pit<-ptt<-rep(0,20)
+     pct[1]<-p[1]
+     pit[1]<-p[2]
+     ptt[1]<-1-p[1]-p[2]
+    for(i in 2:20){
+     pc.old<-pct[i-1]
+     pi.old<-pit[i-1]
+     pt.old<-ptt[i-1]
+         
+     den<-pc.old^2+2*pc.old*pi.old+2*pc.old*pt.old
+     ncc<-nc*pc.old^2/den
+     nci<-2*nc*pc.old*pi.old/den
+     nct<-2*nc*pc.old*pt.old/den
+     nii<-ni*pi.old^2/(pi.old^2+2*pi.old*pt.old)
+     nit<-2*ni*pi.old*pt.old/(pi.old^2+2*pi.old*pt.old)
+      
+     pct<-(2*ncc+nci+nct)/(2*n)
+     pit<-(2*nii+nit+nci)/(2*n)
+     ptt<-(2*ntt+nct+nit)/(2*n)
+     }
+     return(list(pct=pct,pit=pit,ptt=ptt))
+     }
>  
>    n.obs<-c(85,196,341) # observed data,n_c,n_I,n_T
>    p<-c(1/3,1/3)
>    a<-moth(p,n.obs)
0.3333333 0.3333333
>    pct<-a$pct
>    pit<-a$pit
>    ptt<-a$ptt
>    #convergence diagnostics
>    # statistic R
>    rcc=sqrt( (diff(pct)^2+diff(pit)^2)/(pct[-20]^2+pit[-20]^2) )
>    rcc=c(0,rcc)  #adjusts the length to make the table below
>      
>    d1=(pct[-1]-pct[20])/(pct[-20]-pct[20])
>    d1=c(d1,0)
>    d2=(pit[-1]-pit[20])/(pit[-20]-pit[20])
>    d2=c(d2,0)
>
>    #Table output
>    print(cbind(pct,pit,rcc,d1,d2)[1:9,],digits=5)
           pct     pit        rcc       d1      d2
[1,] 0.333333 0.33333 0.0000e+00 0.042502 0.33659
[2,] 0.081994 0.23741 5.7069e-01 0.036933 0.18765
[3,] 0.071249 0.19787 1.6312e-01 0.036727 0.17780
[4,] 0.070852 0.19036 3.5756e-02 0.036719 0.17624
[5,] 0.070837 0.18902 6.5860e-03 0.036719 0.17595
[6,] 0.070837 0.18879 1.1683e-03 0.036719 0.17589
[7,] 0.070837 0.18875 2.0580e-04 0.036719 0.17588
[8,] 0.070837 0.18874 3.6205e-05 0.036720 0.17587
[9,] 0.070837 0.18874 6.3679e-06 0.036736 0.17587
>
> ### Example 2 Bayes EM: mode of posterior
> #generate observed data
> theta=0.3
> n<-50
> x<-drop(rmultinom(1,n,c((2+theta)/4,(1-theta)/2,theta/4)))
> # prior paramters
> a<-.01
> b<-.01
>
>
> th<-rep(0,20)
> th[1]<-0.2
> for(t in 2:20){
+   num<-x[1]*th[t-1]/(2+th[t-1])+x[3]+a-1
+   den<-x[1]*th[t-1]/(2+th[t-1])+x[2]+x[3]+a+b-2
+   th[t]<-num/den
+   }
> rcc<-sqrt( (diff(th)^2+diff(th)^2)/(th[-20]^2+th[-20]^2) )  
> print(cbind(th,c(0,rcc)),digits=5)
           th           
[1,] 0.20000 0.0000e+00
[2,] 0.19297 3.5148e-02
[3,] 0.18928 1.9131e-02
[4,] 0.18732 1.0366e-02
[5,] 0.18627 5.6016e-03
[6,] 0.18570 3.0227e-03
[7,] 0.18540 1.6298e-03
[8,] 0.18524 8.7843e-04
[9,] 0.18515 4.7333e-04
[10,] 0.18510 2.5502e-04
[11,] 0.18508 1.3739e-04
[12,] 0.18506 7.4013e-05
[13,] 0.18506 3.9871e-05
[14,] 0.18505 2.1478e-05
[15,] 0.18505 1.1570e-05
[16,] 0.18505 6.2329e-06
[17,] 0.18505 3.3576e-06
[18,] 0.18505 1.8087e-06
[19,] 0.18505 9.7435e-07
[20,] 0.18505 5.2488e-07
> ### Example 3 SEM algorithm
>
> #First, run EM using the code in the basic EM algorithm given in Example 1.
> #Then we run SEM as follows.  
>
> #Starting values for SEM
> #Choosing the parameter values at the 20th iteration of EM.  
> pcthat=pct[20]   
> pithat=pit[20]
> ptthat=ptt[20]
>
> #Initialize the parameter estimate vectors, one vector for each phenotype
> sempct=rep(0,20)  
> sempit=rep(0,20)  
> semptt=rep(0,20)
> sempct[1]=0.07      #Initial SEM estimates of the parameters
> sempit[1]=0.19
> semptt[1]=0.74
> rij=array(0,c(3,3,20))
>
> #This for loop implements the SEM algorithm for the peppered moth
> #example.  
>
> #Below t is the number of SEM iterations
> for (t in 2:20) {
+   #Take standard EM step
+   denom1=(sempct[t-1]^2+2*sempct[t-1]*sempit[t-1]+2*sempct[t-1]*semptt[t-1])
+   ncct=nc*sempct[t-1]^2/denom1
+   ncit=2*nc*sempct[t-1]*sempit[t-1]/denom1
+   nctt=nc-ncct-ncit
+   niit=ni*sempit[t-1]^2/(sempit[t-1]^2+2*sempit[t-1]*semptt[t-1])
+   nitt=ni-niit
+   nttt=nt
+   sempct[t]=(2*ncct+ncit+nctt)/(2*n)
+   sempit[t]=(2*niit+ncit+nitt)/(2*n)
+   semptt[t]=(2*nttt+nctt+nitt)/(2*n)
+   #SEM
+   #loop over each parameter
+   for (j in 1:3) {
+     #start at estimates from the the 20th iteration of EM
+     sempj=c(pcthat,pithat,ptthat)  
+
+     #replace the jth element of sempj with the most recent EM estimate
+     sempj[j]=c(sempct[t],sempit[t],semptt[t])[j]
+
+     #Take one EM step for sempj   
+     denom1=(sempj[1]^2+2*sempj[1]*sempj[2]+2*sempj[1]*sempj[3])
+     ncct=nc*sempj[1]^2/denom1
+     ncit=2*nc*sempj[1]*sempj[2]/denom1
+     nctt=nc-ncct-ncit
+     niit=ni*sempj[2]^2/(sempj[2]^2+2*sempj[2]*sempj[3])
+     nitt=ni-niit
+     nttt=nt
+     nextstep=c((2*ncct+ncit+nctt)/(2*n),(2*niit+ncit+nitt)/(2*n),
+        (2*nttt+nctt+nitt)/(2*n))
+
+     # Calculate rij.  This is (4.44) on page 101
+     rij[,j,t]=(nextstep-c(pcthat,pithat,ptthat))/
+                   (sempj[j]-c(pcthat,pithat,ptthat)[j])
+       }
+   }
错误: 找不到对象'nc'
>   for(t in 1:20){
+    cat(t,sempct[t],sempit[t],semptt[t])
+    cat("\n")
+    print(rij[,,t])
+    cat("\n")
+    }
1 0.07 0.19 0.74
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

2 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

3 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

4 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

5 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

6 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

7 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

8 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

9 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

10 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

11 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

12 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

13 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

14 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

15 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

16 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

17 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

18 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

19 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

20 0 0 0
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0

> #Note that the algorithm becomes unstable on 8th iteration
> #Below is the output for iteration 7
> #EM after 20 iterations
> #> cbind(pct,pit,ptt)[20,]
> #       pct        pit        ptt
> #0.07083691 0.18873652 0.74042657
>
> #SEM after 7 iterations (the way the indices work, this is index 6
> #> cbind(sempct,sempit,semptt)[6,]
> #    sempct     sempit     semptt
> #0.07083691 0.18873670 0.74042639
>
> rij[,,6]
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    0
[3,]    0    0    0
> #             [,1]         [,2]        [,3]
> #[1,]  0.034117920 -0.002601059 -0.00260106
> #[2,] -0.006910223  0.140676982 -0.03519589
> #[3,] -0.027207449 -0.138075923  0.03779695
>
> #Now need iyhat (inverse of the information matrix for Y)  used in (4.43)
>
> #Equations (4.5)-(4.9) on page 93
> denom1=(pcthat^2+2*pcthat*pithat+2*pcthat*ptthat)
> ncct=nc*pcthat^2/denom1
错误: 找不到对象'nc'
> ncit=2*nc*pcthat*pithat/denom1
错误: 找不到对象'nc'
> nctt=nc-ncct-ncit
错误: 找不到对象'nc'
> niit=ni*pithat^2/(pithat^2+2*pithat*ptthat)
错误: 找不到对象'ni'
> nitt=ni-niit
错误: 找不到对象'ni'
> nttt=nt
错误: 找不到对象'nt'
>

> #Required derivatives for iyhat (used in eqn (4.41) on page 101)
> d20q=-(2*ncct+ncit+nctt)/pcthat^2 -(2*nttt+nitt+nctt)/(ptthat^2)
错误: 找不到对象'ncct'
> d02q=-(2*niit+ncit+nitt)/pithat^2 -(2*nttt+nitt+nctt)/(ptthat^2)
错误: 找不到对象'niit'
> d12q=-(2*nttt+nitt+nctt)/(ptthat^2)
错误: 找不到对象'nttt'
>  
> iyhat=-cbind(c(d20q,d12q),c(d12q,d02q))
错误于cbind(c(d20q, d12q), c(d12q, d02q)) : 找不到对象'd20q'
>
> solve(iyhat)
错误于solve(iyhat) : 找不到对象'iyhat'
> #[1,]  5.290920e-05 -1.074720e-05
> #[2,] -1.074720e-05  1.230828e-04
>
> #Since the percentages of the 3 moth phenotypes add up to 1, we only
> #presented results for carbonaria and insularia phenotypes.
>
> psiprime=rij[,,6] #(7th iteration)
> psiprime22=psiprime[-3,-3]
>
> #variance matarix of MLE
> varhat=solve(iyhat)%*%solve(diag(2)-t(psiprime22))
错误于solve(iyhat) : 找不到对象'iyhat'
> varhat=(varhat+t(varhat))/2
错误: 找不到对象'varhat'
> varhat
错误: 找不到对象'varhat'
> #              [,1]          [,2]
> #[1,]  5.481298e-05 -1.223007e-05
> #[2,] -1.223007e-05  1.433249e-04
>
> #sd(pc), sd(pi)
> sqrt(diag(varhat))              
错误于diag(varhat) : 找不到对象'varhat'
> #[1] 0.007403579 0.011971838
>
>
> #cor(pc,pi)
> varhat[1,2]/prod(sqrt(diag(varhat)))  
错误: 找不到对象'varhat'
> #[1] -0.1379833
>
> #var(pt)
> varhat[1,1]+varhat[2,2]+2*varhat[1,2]  
错误: 找不到对象'varhat'
> #[1] 0.0001736777
>
> #sd(pt)
> sqrt(sum(varhat))
错误: 找不到对象'varhat'
> #[1] 0.01317868
>
> #cor(pc,pt)
> (-varhat[1,1]-varhat[1,2])/(sqrt(varhat[1,1])*sqrt(sum(varhat)))
错误: 找不到对象'varhat'
> #[1] -0.4364370
>
> #cor(pi,pt)
> (-varhat[2,2]-varhat[1,2])/(sqrt(varhat[2,2])*sqrt(sum(varhat)))
错误: 找不到对象'varhat'
> #[1] -0.8309075
>

epoh老师，您好！
        您看一下附件中的EGARCH程序是不是用jplv7中的ucsd-garch工具箱可以运行？
我运行总是出错，您看如何运行？
非常感谢！

怎么也没找到Gauss的code,但是找到了一个Matlab code.是一篇德国人写的论文，现介绍了模型，然后还应用了德国的数据，最后提供了Matlab code。大家看看吧~真遗憾，Gauss 竟然没有现成的程序。
Markov-switching EGARCH (Matlab).pdf

本文来自: 人大经济论坛 Gauss专版 版，详细出处参考： https://bbs.pinggu.org/forum.php? ... 1&from^^uid=11232

z兄，fmincon是matlab自带的
从help也可以看出，如果不是自带的，是没有详细的例子的

上面EGARCH里面有很多m文件，你分开运行过吗？如果没有，那是需要分开的

EM.R

少定义nc,ni.nt

请自行补上.

另两个问题,ywh兄所言甚是.

### Example 3 SEM algorithm

#First, run EM using the code in the basic EM algorithm given in Example 1.

#Then we run SEM as follows.  

n.obs<-c(85,196,341)

nc<-n.obs[1]

ni<-n.obs[2]

nt<-n.obs[3]

#Starting values for SEM

epoh老师，您好！
    附件中的fortran程序如何运行？需要下载什么软件才能运行firtran程序？
另外，附件中的matlab程序运行出现以下错误？为什么？如何修改？

1\Simulatione.m

simulation

s =

     1

??? Undefined function or variable "np".

Error in ==> simulatione at 41

year=zeros(np*TS,1);

>>

2\solvee.m

Warning: Input arguments must be scalar.

> In transzx at 3

  In solvee at 8

Warning: Size vector should be a row vector with integer elements.

> In transzx at 3

  In solvee at 8

Warning: Input arguments must be scalar.

> In transzx at 3

  In solvee at 8

Warning: Size vector should be a row vector with integer elements.

> In transzx at 3

  In solvee at 8

??? Undefined function or variable "Tx11".

Error in ==> transzx at 10

Tzx11=Tz.*Tx11;

Error in ==> solvee at 8

transzx;

>>

3\xprobit.m

??? Undefined function or variable "nkgrid".

Error in ==> xprofit at 5

pix=zeros(ngrid,nkgrid);

>>

Dprobit.m

??? Undefined function or variable "nkgrid".

Error in ==> dprofit at 6

pid=zeros(ngrid,nkgrid);

>>

Parameters.m

运行没有任何反应？

还有cd.doc中的stata错误如何修改？
非常感谢！

很久没用fortran了

由于mcmc.f90要用到IMSL.mod,numerical_libraries.mod

PROGRAM mcmc

!=======================

!Use these modules

!=======================

USE IMSL

USE numerical_libraries

USE parameterse

单独抓IMSL.mod,numerical_libraries.mod都失败

最后使用内含IMSL Fortran Library的

Compaq Visual Fortran 6.6 Pro

才顺利Build mcmc.exe

mcmc_fortran.rar

  

double click mcmc.exe

执行时会自动产生2个文件

  accept.txt,

  param.txt,

zhangtao兄

Rdexp_summary内的readme

已经说明执行function是main.m

依序执行load_data,parameters,.....

如此才不会产生错误.

%LOAD DATA

load_data;

%DEFINE PARAMETERS GIVEN

parameters;

%DEFINE GRIDS OF CAPITAL AND PRODUCTIVITY (X)

gridkx;

......

......

epoh老师，您好！
       我想问一下：如何把R中的tsDyn这个包载入到splus中进行运用和进行相关的
计量分析？

哈哈!这个有困难.

> setwd('C:/R_threg')
> melanoma<-read.table("melanoma.txt",header=TRUE)
> threg(data=melanoma,lny0=age,mu=thick+age+sex,timevar=survtime,failure=status)
错误: 没有"threg"这个函数
>

epoh老师，您好！
    附件中R程序运行出现以上问题，您看如何解决？
非常感谢！

另外，tsDyn包载入splus中运行，您说有困难？那困难在什么地方？
能实现吗？如果能实现，如何实现？

melanoma<-read.table("melanoma.txt",header=TRUE)
source("threg.R")
threg(data=melanoma,lny0=age,mu=thick+age+sex,timevar=survtime,failure=status)
               Coef. Est. Std. Err.      z  P>|z|
lny0: age       0.0033790 0.0040379  0.837 0.4026
lny0: constant  3.5621735 0.2288590 15.565 0.0000
mu: sex        -0.0079129 0.0045327 -1.746 0.0808
mu: age        -0.0003492 0.0002030 -1.720 0.0854
mu: thick      -0.0029611 0.0007753 -3.819 0.0001
mu: constant    0.0370407 0.0105242  3.520 0.0004
#######

Q2:tsDyn包载入splus中运行，有困难？困难在什么地方?

   你是否还记得package"pkgutils"

   s-plus就是用pkgutils,

   建构供s-plus使用的package

   目前R与S-plus语法已经相差很多

   要用pkgutils建构符合s-plus用的

   tsDyn package 肯定会出错.

遇到同样的问题，正寻找用lingo做DEA的教程呢？有没有大虾分享或推荐一下呀。。

QIANG`````````

epoh老师好，我想问问你一个关于DYNARE问题，我编好程序后，开始运行，出现了错误提示：
Starting Dynare (version 4.2.1).
Starting preprocessing of the model file ...
ERROR: fs2000.mod:96.71: syntax error, unexpected ';'

??? Error using ==> dynare
DYNARE: preprocessing failed


这是我编好的程序，不知道是哪里出错了：
var c l y k i gi gc z b rk v r f kg a sb w taul tauk tauc ub ul ua uz utl utk utc ugc;
varexo e_b e_l e_a1 e_a2 e_z e_tl e_tk e_tc e_gc;

parameters alp bet gam h kappa phi0 phi1 phi2 phi3 phi4 phi5 phi6 phi7 phi8 phi9 phi10 phi11 del0 del1 del2 delg alpg rhob rhol rhoa1 rhoa2 rhoz rhotl rhotk rhotc rhogc sigmab sigmal sigmaa1 sigmaa2 sigmaz sigmatk sigmatl sigmatc sigmagc psiz psitk psitl gammaz gammatk gammatl gammagc;

alp = 0.436;
bet = 0.9897;
gam = 0.87;
h = 0.7;
kappa = 1.1;
phi0 = 0;
phi1 = 0.03125;
phi2 = 0.03125;
phi3 = 0.03125;
phi4 = 0.09375;
phi5 = 0.09375;
phi6 = 0.09375;
phi7 = 0.09375;
phi8 = 0.09375;
phi9 = 0.09375;
phi10 = 0.09375;
phi11 = 0.09375;
del0 = 0.0274;
del1 = 0.035;
del2 = 0.7;
delg = 0.025;
alpg = 0.2;
rhob = 0.5;
rhol = 0.5;
rhoa1 = 0.5;
rhoa2 = 0.5;
rhoz = 0.5;
rhotl = 0.5;
rhotk = 0.5;
rhotc = 0.5;
rhogc = 0.5;
sigmab = 1;
sigmal = 1;
sigmaa1 = 1;
sigmaa2 = 1;
sigmaz = 1;
sigmatl = 1;
sigmatk = 1;
sigmatc = 1;
sigmagc = 1;
psiz = 0.2;
psitk = 1;
psitl = 0.5;
gammaz = 0.15;
gammatk = 0.15;
gammatl = 0.15;
gammagc = 0.15;




model;
(c-h*c)^(-gam) = ((1+tauc)/ub)*bet*r*ub(+1)*(c(+1)-h*c)^(-gam)/(1+tauc(+1));
l^(1+kappa) = (1/(ul*(1+tauc))*((c-h*c(-1))^(-gam))*(1-taul)*(1-alp)*y;
y = rk*v*k(-1)+w*l;
ub = rhob*ub(-1)+sigmab*e_b;
ul = rhol*ul(-1)+sigmal*e_l;
k = (1-f)*k(-1)+i;
f = del0+del1*(v-1)+(del2*(v-1)^2)/2;
y = ua*((v*k(-1))^alp)*(i^(1-alp))*kg(-1)^alpg;
ua = rhoa1*ua(-1)+sigmaa1*e_a1;
i = del0*k;
z = tauc*c+tauk*rk*v*k(-1)+taul*w*l+b-gc-gi-r(-1)*b(-1);
rk = alp*y/k(-1);
w = (1-alp)*/l;
kg(-1) = (1-delg)*kg(-2)+a(-12);
gi = phi0*a+phi1*a(-1)+phi2*a(-2)+phi3*a(-3)+phi4*a(-4)+phi5*a(-5)+phi6*a(-6)+phi7*a(-7)+phi8*a(-8)+phi9*a(-9)+phi10*a(-10)+phi11*a(-11);
a = rhoa2*a(-1)+sigmaa2*e_a2;
y = c+i+gi+gc;
zhat = -psia*yhat-gammaz*sbhat(-8)+uzhat;
uzhat = rhoz*uzhat(-1)+sigmaz*e_z;
taulhat = psitl*yhat+gammatl*sbhat(-8)+utlhat;
utlhat = rhotl*utlhat(-1)+sigmatl*e_tl;
taukhat = psitk*yhat+gammatk*sbhat(-8)+utkhat;
utkhat = rhotk*utkhat(-1)+sigmatk*e_tk;
tauchat = utchat;
utchat = rhotc*utchat(-1)+sigmatc*e_tc;
gchat = -gammagc*sbhat(-8)+ugchat;
ugchat = rhogc*ugchat(-1)+sigmagc*e_gc;
sb = b/y;
end;

initval;
sb = 0.04;
v = 1;
rk = 0.057;
i =0.2096;
c = 0.45;
k = 7.649;
l = 0.2075;
y = 1;
gc = 0.0713;
b = 0.04;
gi = 0.0028;
kg = 3.627;
z = 0.03949
w = 2.7187;
taul = 0.1;
tauk = 0.1;
tauc = 0.1;
uz = 0;
utk = 0;
utl = 0;
utc = tauc;
ugc = 0;
ub = 1;
ua = 1.7084;
ul = 0.03674;
a = gi;
f = del0;
r = 1/bet;
e_b = 0;
e_l = 0;
e_a1 = 0;
e_a2 = 0;
e_z = 0;
e_tl = 0;
e_tk = 0;
e_tc = 0;
e_gc = 0;
end;


shocks;
var e_b; stderr 0.014;
var e_l; stderr 0.005;
var e_a1; stderr 0.005;
var e_a2; stderr 0.005;
var e_z; stderr 0.005;
var e_tl; stderr 0.005;
var e_tk; stderr 0.005;
var e_tc; stderr 0.005;
var e_gc; stderr 0.005;
end;

steady(solve_algo=0);

check;

estimated_params;
alp, beta_pdf, 0.4, 0.02;
bet, beta_pdf, 0.99, 0.002;
gam, gamma_pdf, 0.87, 0.5;
kappa , gamma_pdf, 1.1, 0.5;
h, beta_pdf, 0.7, 0.2;
del2, beta_pdf, 0.7, 0.5;
rhob, beta_pdf, 0.5, 0.2;
rhol, beta_pdf, 0.5, 0.2;
rhoa1, beta_pdf, 0.5, 0.2;
rhoa2, beta_pdf, 0.5, 0.2;
rhoz, beta_pdf, 0.5, 0.2;
rhotl, beta_pdf, 0.5, 0.2;
rhotk, beta_pdf, 0.5, 0.2;
rhotc, beta_pdf, 0.5, 0.2;
rhogc, beta_pdf, 0.5, 0.2;
sigmab, inv_gamma_pdf, 1, inf;
sigmal, inv_gamma_pdf, 1, inf;
sigmaa1, inv_gamma_pdf, 1, inf;
sigmaa2, inv_gamma_pdf, 1, inf;
sigmaz, inv_gamma_pdf, 1, inf;
sigmatl, inv_gamma_pdf, 1, inf;
sigmatk, inv_gamma_pdf, 1, inf;
sigmatc, inv_gamma_pdf, 1, inf;
sigmagc, inv_gamma_pdf, 1, inf;
psiz, gamma_pdf, 0.2, 0.1;
psitk, gamma_pdf, 1, 0.3;
psitl, gamma_pdf, 0.5, 0.25;
gammaz, normal_pdf, 0.15, 0.1;
gammatk, normal_pdf, 0.15, 0.1;
gammatl, normal_pdf, 0.15, 0.1;
gammagc, normal_pdf, 0.15, 0.1;
stderr e_b, inv_gamma_pdf, 0.035449, inf;
stderr e_l, inv_gamma_pdf, 0.008862, inf;
stderr e_a1, inv_gamma_pdf, 0.008862, inf;
stderr e_a2, inv_gamma_pdf, 0.008862, inf;
stderr e_z, inv_gamma_pdf, 0.008862, inf;
stderr e_tl, inv_gamma_pdf, 0.008862, inf;
stderr e_tk, inv_gamma_pdf, 0.008862, inf;
stderr e_tc, inv_gamma_pdf, 0.008862, inf;
stderr e_gc, inv_gamma_pdf, 0.008862, inf;
end;

varobs y c l;

estimation(order=1, datafile=fsdat_simul, mh_replic=2000, mh_nblocks=2, mh_jscale=0.8);


/*
* The following lines were used to generate the data file. If you want to
* generate another random data file, comment the "estimation" line and uncomment
* the following lines.
*/

//stoch_simul(periods=200, order=1);
//datatomfile('fsdat_simul', char('y', 'c', 'l'));

谢谢epoh老师!

model部份,

初看有下列3点错

可能还有,请再仔细检查.

1.

l^(1+kappa) = (1/(ul*(1+tauc))*((c-h*c(-1))^(-gam))*(1-taul)*(1-alp)*y;

  漏掉 " ) "

我先改成这样:

l^(1+kappa) = (1/(ul*(1+tauc)))*((c-h*c(-1))^(-gam))*(1-taul)*(1-alp)*y;

%%%%%%

2.

w = (1-alp)* ? /l;  

%%%%%%

3.

zhat : Unknown symbol

谢谢epoh老师，的确是这里出了错误！修改后，就好了。但是错误一直像冒水泡一样的不断出现。。。无语了，还得请教你这个问题：

Starting Dynare (version 4.2.1).
Starting preprocessing of the model file ...
Substitution of endo lags >= 2: added 19 auxiliary variables and equations.
Found 47 equation(s).
Evaluating expressions...done
Computing static model derivatives:
- order 1
Computing dynamic model derivatives:
- order 1
Processing outputs ...done
Preprocessing completed.
Starting MATLAB/Octave computing.

Residuals of the static equations:

Equation number 1 : 0
Equation number 2 : -1.2211
Equation number 3 : 0.15673
Equation number 4 : 0.5
Equation number 5 : 0.5
Equation number 6 : 0
Equation number 7 : 0.06635
Equation number 8 : 0.49994
Equation number 9 : 0.5
Equation number 10 : 0.033175
Equation number 11 : 0.065673
Equation number 12 : -0.31346
Equation number 13 : 0
Equation number 14 : 0
Equation number 15 : 0
Equation number 16 : 0
Equation number 17 : -0.09375
Equation number 18 : -0.00114
Equation number 19 : 0
Equation number 20 : 0.05896
Equation number 21 : 0
Equation number 22 : 0.06886
Equation number 23 : 0
Equation number 24 : -0.05
Equation number 25 : 0.075
Equation number 26 : 0
Equation number 27 : 0
Equation number 28 : 0


??? Error using ==> steady_
The steadystate values returned by fs2000_steadystate.m don't solve the static model!

Error in ==> steady at 54
steady_;

Error in ==> fs2000 at 573
steady;

Error in ==> dynare at 132
evalin('base',fname) ;

下面是设置的steady state的程序：
function [ys,check] = fs2000_steadystate(ys,exe)
  global M_
  
  alp = M_.params(1);
  bet = M_.params(2);
  gam = M_.params(3);
  h = M_.params(4);
  kappa = M_.params(5);
  phi0 = M_.params(6);
  phi1 = M_.params(7);
  phi2 = M_.params(8);
  phi3 = M_.params(9);
  phi4 = M_.params(10);
  phi5 = M_.params(11);
  phi6 = M_.params(12);
  phi7 = M_.params(13);
  phi8 = M_.params(14);
  phi9 = M_.params(15);
  phi10 = M_.params(16);
  phi11 = M_.params(17);
  del0 = M_.params(18);
  del1 = M_.params(19);
  del2 = M_.params(20);
  delg = M_.params(21);
  alpg = M_.params(22);
  rhob = M_.params(23);
  rhol = M_.params(24);
  rhoa1 = M_.params(25);
  rhoa2 = M_.params(26);
  rhoz = M_.params(27);
  rhotl = M_.params(28);
  rhotk = M_.params(29);
  rhotc = M_.params(30);
  rhogc = M_.params(31);
  sigmab = M_.params(32);
  sigmal = M_.params(33);
  sigmaa1 = M_.params(34);
  sigmaa2 = M_.params(35);
  sigmaz = M_.params(36);
  sigmatl = M_.params(37);
  sigmatk = M_.params(38);
  sigmatc = M_.params(39);
  sigmagc = M_.params(40);
  psiz = M_.params(41);
  psitk = M_.params(42);
  psitl = M_.params(43);
  gammaz = M_.params(44);
  gammatk = M_.params(45);
  gammatl = M_.params(46);
  gammagc = M_.params(47);
  taul = M_.params(48);
  tauk = M_.params(49);
  tauc = M_.params(50);
  check = 0;
  
  
  
  
  v = 1;
  r = 1/bet;
  rk = (1/bet-(1-del0))/(1-tauk);
  y = 0.5;
  k = 0.5;
  c = 0.5;
  l = 0.5;
  sb = 0.3396;
  ub = 1;
  ul = 1;
  ua = 1;
  uz = 0;
  utl = 0;
  ugc = 0;
  utk = 0;
  
  
  i = del0*k;
  w = (1-alp)*y/l;
  kg =(del0^(1-alp)*k/y)^(1/alpg);
  gi = delg*kg;
  a = gi;
  gc= c+i+gi-y;
  b = sb*y;   
  z = tauc*c+tauk*alp*y+taul*w*l-gi-gc+(1-r)*b;
  f = del0;
  yhat = (taul-gammatl*sb)/psitl;
  zhat = -psiz*yhat-gammaz*sb;
  gchat = -gammagc*sb;
  utc = tauc;
  
  
  
  
  
  ys =[
c
l
y
k
i
gi
gc
z
b
rk
v
r
f
kg
a
sb
w
zhat
yhat
gchat
ub
ul
ua
uz
utl
utk
utc
ugc    ];
不知道steady 的程序哪里出了问题，还请你帮忙给检查一下，谢谢你！

哈哈!＂美＂的冒泡

问题出在这里,有点麻烦.

方便的话你把修改过的

model代码,放在我的短信息.

(有变更的部分)

我来实际跑程序,

这样好跟你的结果一致.

%%%%%%%

要不你就跟底下这篇文献一样

不要用fs2000_steadystate试试

Stochastic simulation of a DSGE model for Brazil.pdf
Brazilmodel_v4.mod
data33.mat

   

epoh 老师呀，你竟然说美的泡泡，真是太佩服你了！ 我是心里堵得慌呀，不敢面对那些“泡泡”了。。。
我按照您给的文献修改了一下，不用steadystate了 ,结果还是不断冒泡呀，实在是快被淹窒息了。
现将这些代码发给您，请您救我呀！太感激您了！！！谢谢！

Residuals of the static equations:

Equation number 1 : 0
Equation number 2 : -1.2211
Equation number 3 : 0.15673
Equation number 4 : 0
Equation number 5 : 0.06635
Equation number 6 : 0.49994
Equation number 7 : 0.033175
Equation number 8 : 0.065673
Equation number 9 : -0.31346
Equation number 10 : 0
Equation number 11 : 0
Equation number 12 : 0
Equation number 13 : 0
Equation number 14 : -0.09375
Equation number 15 : 0
Equation number 16 : -0.00114
Equation number 17 : 0.05896
Equation number 18 : 0.06886
Equation number 19 : -0.05
Equation number 20 : 0
Equation number 21 : 0.5
Equation number 22 : 0.5
Equation number 23 : 0.5
Equation number 24 : 0
Equation number 25 : 0
Equation number 26 : 0
Equation number 27 : 0.075
Equation number 28 : 0
************

solve_algo = 0

solve_algo = 1

solve_algo = 2

solve_algo = 3

我都试了,Residuals都很大

要麻烦你检查一下

1.你手算的 initial values 对不对?

2.每条输入的equation对不对?

**********

governinv_dynamic.m
% Model equations
%

residual = zeros(47, 1);
.....
.....
lhs =T15;
rhs =T25*T29/(1+params(8));
residual(1)= lhs-rhs;
lhs =y(36)^(1+params(5));
rhs =T48*y(37);
residual(2)= lhs-rhs;
lhs =y(37);
rhs =y(44)*y(45)*y(2)+y(36)*y(50);
residual(3)= lhs-rhs;
.....
.....
lhs =y(79);
rhs =y(32);
residual(46)= lhs-rhs;
lhs =y(80);
rhs =y(33);
residual(47)= lhs-rhs;

设定了几个初始值，然后其它的 initial values就是手算的。

Jpfeifer认为：
“ If you computed the theoretical steady state, put the formulas in the initval-block and put resid(1) before the steady command. If the residuals are not 0, your theoretical steady state is wrong.”不是很明白，epoh老师能理解是什么意思吗？
错误就集中在我输入的equations 和initial values的正确性了。
谢谢epoh老师对我的帮助！ you are so kind! thank you !

put resid(1) before the steady command.

就是:

initval;
v = 1;
rk = 0.042;
.....
.....
e_tk = 0;
e_tc = 0;
e_gc = 0;
end;

resid(1);
/*
steady;

check;

estimated_params;
alp, beta_pdf, 0.4, 0.02;
.....
.....
estimation(datafile=gidat_simul,.....);

*/
%%%%%%%%%%%%%%%%

dynare governinv
Substitution of endo lags >= 2: added 19 auxiliary variables and equations.???
Found 47 equation(s).

Residuals of the static equations:

Equation number 1 : 0
Equation number 2 : -1.2211
Equation number 3 : 0.15673
Equation number 4 : 0
Equation number 5 : 0.06635
Equation number 6 : 0.49994
Equation number 7 : 0.033175
Equation number 8 : 0.065673
Equation number 9 : -0.31346
Equation number 10 : 0
Equation number 11 : 0
Equation number 12 : 0
Equation number 13 : 0
Equation number 14 : -0.09375
Equation number 15 : 0
Equation number 16 : -0.00114
Equation number 17 : 0.05896
Equation number 18 : 0.06886
Equation number 19 : -0.05
Equation number 20 : 0
Equation number 21 : 0.5
Equation number 22 : 0.5
Equation number 23 : 0.5
Equation number 24 : 0
Equation number 25 : 0
Equation number 26 : 0
Equation number 27 : 0.075
Equation number 28 : 0

%%%%%%%%%%%
dynare Brazilmodel_v4

Starting Dynare (version 4.2.2).
Starting preprocessing of the model file ...
Found 25 equation(s).
Evaluating expressions...done
Computing static model derivatives:
- order 1
Computing dynamic model derivatives:
- order 1
- order 2
Processing outputs ...done
Preprocessing completed.
Starting MATLAB/Octave computing.

Residuals of the static equations:

Equation number 1 : 0
Equation number 2 : 0
Equation number 3 : 0
Equation number 4 : 0
Equation number 5 : 0
Equation number 6 : 0
Equation number 7 : 0
Equation number 8 : 0
Equation number 9 : 0
Equation number 10 : 0
Equation number 11 : 0
Equation number 12 : 0
Equation number 13 : 0
Equation number 14 : 0
Equation number 15 : 0
Equation number 16 : 0
Equation number 17 : 0
Equation number 18 : 0
Equation number 19 : 0
Equation number 20 : 0
Equation number 21 : 0
Equation number 22 : 0
Equation number 23 : 0
Equation number 24 : 0
Equation number 25 : 0

%%%%%%%%%

你的是啥模型?

Dynare v4 - User Guide
page 35/106
For complicated models,finding suitable initial values for
the endogenous variables is the trickiest part of finding the
equilibrium of that model. Often,it is better to
start with a smaller model and add new variables one by one.

epoh老师太精通了，我是想把这两篇文章中的model做成DYNARE的程序，自己编的，是不是有很大的问题，我不知道选择model的标准是什么，我是以有多少内生变量有多少方程来编的，哎，我好像是什么都不懂呢~~~~