adp.zip_ADP算法_ADP近似动态规划在经济学中的应用_adp

clear all; close all; % Define the structural parameters of the model, % i.e. policy invariant preference and technology parameters % alpha : capital's share of output % beta : time discount factor % delta : depreciation rate 折旧率 % sigma : risk-aversion parameter, also intertemp. subst. param. % gamma : unconditional expectation of the technology parameter alpha = .35; beta = .98; delta = .025; sigma = 2; zbar = 5; % Find the steady-state level of capital as a function of % the structural parameters 稳态水平作为结构参数的函数 kstar = ((1/beta - 1 + delta)/(alpha*zbar))^(1/(alpha-1)); % Define the number of discrete values k can take gk = 101; k = linspace(0.95*kstar,1.05*kstar,gk); % Compute a (gk x gk) dimensional consumption matrix c % for all the (gk x gk) values of k_t, k_t+1 for h = 1 : gk for i = 1 : gk c(h,i) = zbar*k(h)^alpha + (1-delta)*k(h) - k(i); if c(h,i) < 0 c(h,i) = 0; end % h is the counter for the endogenous state variable k_t % i is the counter for the control variable k_t+1 end end % Compute a (gk x gk) dimensional consumption matrix u % for all the (gk x gk) values of k_t, k_t+1 for h = 1 : gk for i = 1 : gk if sigma == 1 u(h,i) = log(c(h,i)); else u(h,i) = (c(h,i)^(1-sigma) - 1)/(1-sigma); end % h is the counter for the endogenous state variable k_t % i is the counter for the control variable k_t+1 end end % Define the initial matrix v as a matrix of zeros (could be anything) v = zeros(gk,1); % Set parameters for the loop convcrit = 1E-9; % chosen convergence criterion diff = 1; % arbitrary initial value greater than convcrit iter = 0; % iterations counter while diff > convcrit % for each k_t % find the k_t+1 that maximizes the sum of instantenous utility and % discounted continuation utility for h = 1 : gk Tv(h,1) = max( u(h,:) + beta*v'); end iter = iter + 1; diff = norm(Tv - v); v = Tv; end disp(iter) % Find what gridpoint of the k_t+1 vector which gave the optimal value for h = 1 : gk [Tv,gridpoint] = max(u(h,:) + beta*v'); kgridrule(h) = gridpoint; end % Find what element of k_t+1 which gave the optimal value kdecrule = k(kgridrule); figure plot(k,kdecrule,k,k); xlabel('k_t') ylabel('k_{t+1}') print -djpeg kdecrule.jpg % Compute the optimal decision rule for c as a function of the state % variable for h = 1 : gk cdecrule(h) = zbar*k(h)^alpha + (1-delta)*k(h) - kdecrule(h); end figure plot(k,cdecrule) xlabel('k_t') ylabel('c_t') print -djpeg cdecrule.jpg