ECON 474-625: Numerical Methods for Economists
Assignment 3: Unconstrained Optimization Methods
(Due date: Sunday February 25 before 11:30pm in Learn drop box)
This is an individual assignment. You can help each other, but the final result must be yours.
Question 1
In this question, we want to compute the market demand using a two good market and a one dimen- tional unconstrained Newton method. You have to answer the question using expression and analitical derivatives. The market is composed of 1,000 consumers with the following utility function:
To generate you 1000 consumers, you will proceed as folows:
set.seed(seed) s <- runif(1000, 1.5, 5) a1 <- runif(1000, .2, 3) a2 <- runif(1000, 2, 6) I <- runif(1000, 500, 2500) par <- cbind(s=s,a1=a1,a2=a2,I=I)
with the seed number replaced by your student id. Follow the following steps:
a) Write an expression for the utility function and a function that will return either the negative of the utility, the negative of the first derivative or the negative of the second derivarive (with respect to x1) at x1 when prices are p1 and p2, income is I and consumer spends all his income. The function would work as:
getf(ces, list(x1=10, p1=5, p2=7, I=500, a1=1, a2=2, s=3), d=0) ## [1] -200.4395 getf(ces, list(x1=10, p1=5, p2=7, I=500, a1=1, a2=2, s=3), d=1) ## [1] 1.45659 getf(ces, list(x1=10, p1=5, p2=7, I=500, a1=1, a2=2, s=3), d=2) ## [1] 0.04622944
b) Write a function that returns the optimal x1 and x2 using the Newton method. It should solve the following problem:
l <- list(p1=5, p2=7, I=500, a1=1, a2=2, s=3) getsol1(U=ces, x0=3,l=l) ## $sol ## ## 1.734754 70.189461 ## ## $niter ## [1] 7
Here x0 is the starting value for x1 because x2 is replaced by (I −p1x1)/p2, and the optimization is done with respect to x1.
c) The Newton method if not that reliable, expecially when there are restrictions on the possible values for x. For example, if we start too far from the solution, we get the following
getsol1(U=ces, x0=80,l=l,trace=TRUE) ## iteration 1, x1=80, x2=14.29 ## iteration 2, x1=-80.73, x2=129.1 ## Error in if (crit < tol) break: missing value where TRUE/FALSE needed
We can see that the Newton step is too large, which results in x1 being negative in the second step. That could be a problem because we want to solve 1,000 consumer problems with 1,000 different sets of parameters. Choosing a starting value that will make the algorithm converges for all consumers is impossible. In this step, you will modify your function so that it always converge for any sensible starting value (x0 must be less than I/p1). Here is what you will do. If x10 is the initial x1, the Newton method is:
where s = fj(x10)/fjj(x10). If x10 is less than the optimal x1, s will be negative and x10 will increase. We don’t want it to increase too much and make x2 negative. If it happens, we want to lower the value of s. If x10 is greater than the optimal x1, s will be positive and x10 will decrease. If it decreases by too much, it will become negative. We can prevent that by decreasing s. Here is the modification you will add to your function: (i) As long as s > x0, then s = s/2. This will prevent x1 to become negative. (ii) As long as s < (x0 − I/p1), then s = s/2. This will prevent x2 to become negative. Something like the following (before computing the new x1) will do the job:
while (s>x0 | s<(x0-l$I/l$p1)) s <- s/2
The new function should now work for any x0:
getsol2(U=ces, x0=80,l=l) ## $sol ## x1 x2 ## 1.734754 70.189461 ## ## $niter ## [1] 10 getsol2(U=ces, x0=99,l=l) ## $sol ## x1 x2 ## 1.734754 70.189461 ## ## $niter ## [1] 9
Once your function is stable, remove the argument x0 and set the starting value to a sensible value inside the function. Remember that x2 must be positive.
getsol(U=ces, l) ## $sol ## x1 x2 ## 1.734754 70.189461 ## ## $niter ## [1] 8
c) Write a function that computes the total quantity demanded by your market (the sum of all 1,000 consumers) for a given set of prices p = p1, p2 j. Try to make it work for any number of consumers. You would get something like:
library(parallel) marketD(p=c(3,6), par, U=ces, tol=1e-7, maxit=100, mc.cores=4) ## Market demand ## *********************** ## Number of consumers: 1000 ## p1 = 3, p2 = 6 ## x1 = 48968.99, x2 = 222166.76 ## ***********************
Of course, you should all get different values as your par will be different. You are free to create a particular object with its print method if you want, but it is not required.
d)
Write a function that plots the market demand of either x1 for a given p2 or x2 for a given p1. You would get something like:
Hint: If you want to smooth a function with a few points, you can use spline(). For example:
e) Compare the effect of a 20% income on both demands when the other price is 5. I am expecting something like
Question 2
Estimate the ρ from the second assignment using the Newton method and numerical derivative. To do so, you need to have a function that depends on ρ, and returns the SSR. You can use part of your codes from the second assignment. Basically, for each ρ, you need to estimate the model:
set.seed(112233) u <- arima.sim(n=200, model=list(ar=.7)) x1 <- rnorm(200) x2 <- rnorm(200) y <- 1+x1-x2+u X <- cbind(Intercept=1,x1,x2) (ssr <- SSR(rho=.5, y, X, "CO")) ## [1] 202.0671 (ssr <- SSR(rho=.5, y, X, "PW")) ## [1] 202.0674
After, you need a function that computes the first and second derivative of SSR(ρ).
dssr <- dSSR(rho=.5, y, X, "CO",h=.001)) ## [1] -129.2444 (ddssr <- ddSSR(rho=.5, y, X, "PW",h=.001)) ## [1] 709.8033
Then, write the algorithm to minimize SSR. You should get:
getrho(rho0=.1, y, X, "CO") ## $rho ## [1] 0.6819157 ## ## $SSR ## [1] 190.3081 ## ## $niter ## [1] 3 getrho(rho0=.1, y, X, "PW") ## $rho ## [1] 0.6819152 ## ## $SSR ## [1] 190.3086 ## ## $niter ## [1] 3
代写CS&Finance|建模|代码|系统|报告|考试
编程类：C++,JAVA ,数据库,WEB,Linux,Nodejs,JSP,Html,Prolog,Python,Haskell,hadoop算法,系统 机器学习
金融类：统计，计量，风险投资，金融工程，R语言，Python语言，Matlab，建立模型，数据分析，数据处理
服务类：Lab/Assignment/Project/Course/Qzui/Midterm/Final/Exam/Test帮助代写代考辅导
E-mail:850190831@qq.com 微信：BadGeniuscs 工作时间：无休息工作日-早上8点到凌晨3点
如果您用的手机请先保存二维码到手机里面，识别图中二维码。如果用电脑，直接掏出手机果断扫描。
