Applications
% Generation of Tobit Sample rng(2468); nobs = 100; X=[ones(nobs,1),10*rand(nobs,2)]; beta = [5,2,-3]'; epsilon = 3*randn(nobs,1); Y = X * beta + epsilon; % Truncate Ytrunc = (Y > 0) .* Y; %Starting values using OLS beta0 = (X'*X)\(X'*Ytrunc); e0 = Ytrunc-X*beta0; sd = sqrt((e0'* e0)/(nobs-length(beta0))); parm0 = [beta0', sd]; parm0 %Minimization [parm, fval, exitflag, output] = fminsearch(@(parm)tobit_like(parm,Y,X),parm0) %Standard Errors tobit_like2 = @(parm)-tobit_like(parm, Y, X); [hess,err] = hessian(tobit_like2,parm'); varcov = -inv(hess); separm = sqrt(diag(varcov)); tparm = parm' ./ separm; %Printing results fprintf('\nSummary Estimate of Tobit\n') fprintf(' %12s %12s %12s %12s\n',... '','Coef.' , 'Std Error', 'z ') fprintf('%12s %12.4f %12.4f %12.4f\n',... 'const', parm(1,1),separm(1),tparm(1)) fprintf('%12s %12.4f %12.4f %12.4f\n',... 'X2', parm(1,2),separm(2),tparm(2)) fprintf('%12s %12.4f %12.4f %12.4f\n\n',... 'X3', parm(1,3),separm(3),tparm(3)) fprintf('Standard Error of Equation is %7.4f (%7.4f)\n',... parm(1,4),separm(4))
parm0 = 7.4239 1.0862 -1.6571 3.6500 parm = 5.1849 2.0336 -3.1373 3.4256 fval = 145.0262 exitflag = 1 output = iterations: 132 funcCount: 226 algorithm: 'Nelder-Mead simplex direct search' message: 'Optimization terminated:...' Summary Estimate of Tobit Coef. Std Error z const 5.1849 1.2634 4.1038 X2 2.0336 0.1839 11.0551 X3 -3.1373 0.2246 -13.9714 Standard Error of Equation is 3.4256 ( 0.3328)
clc; clear; %Generate data, sample size 20, 100, 500, 2500 rng(56789); x20 = 100*rand(20, 1); x100 = 100*rand(100, 1); x500 = 100*rand(500, 1); x2500 = 100*rand(2500, 1); %Simulation nsimul = 10000; BETA = nan(nsimul, 1); beta = [10; 2] %True value of coefficientsfor s = [20, 100, 500, 2500] switch s %Use different sample size for x case 20 x = x20; case 100 x = x100; case 500 x = x500; case 2500 x = x2500; end X = [ones(size(x,1), 1), x]; % Regressor X for ii = 1:nsimul %simulate nsimul times eps = 20 * randn(size(X,1), 1); %simulated error term y = X * beta + eps; %simulated y value betahat = (X'*X) \ (X'*y); %OLS beta BETA(ii, 1) = betahat(2); %store beta end if s==20 BETA20 = BETA; elseif s==100 BETA100 = BETA; elseif s==500 BETA500 = BETA; elseif s==2500 BETA2500 = BETA; endend fprintf('With 20 obs simulation: mean %6.3f st.dev %6.3f \n',... mean(BETA20), std(BETA20)); fprintf('With 100 obs simulation: mean %6.3f st.dev %6.3f \n',... mean(BETA100), std(BETA100)); fprintf('With 500 obs simulation: mean %6.3f st.dev %6.3f \n',... mean(BETA500), std(BETA500)); fprintf('With 2500 obs simulation: mean %6.3f st.dev %6.3f \n',... mean(BETA2500), std(BETA2500)); n = hist([BETA20,BETA100,BETA500,BETA2500],1.4:0.01:2.6); plot((1.4:0.01:2.6)',n/nsimul,'LineWidth',2); legend('Sample 20','Sample 100','Sample 500','Sample 2500');
beta = 10 2 With 20 obs simulation: mean 2.000 st.dev 0.130 With 100 obs simulation: mean 2.000 st.dev 0.066 With 500 obs simulation: mean 2.001 st.dev 0.032 With 2500 obs simulation: mean 2.000 st.dev 0.014
a = sym('1/3'); b = sym('4/5'); c = sym('c', 'positive'); d = a*b %fraction calculation
d = 4/15
expand( (a-c)^2 ) %expansion
ans = c^2 - (2*c)/3 + 1/9
factor(ans) %factors
ans = [ 1/9, 3*c - 1, 3*c - 1]
syms a b c x y y = a*x^2 + b*x^2 +b*x + c
y = c + b*x + a*x^2 + b*x^2
syms a g g = (3*a)/4 + 3*a/5 - 1.5 + 0.9
g = (27*a)/20 - 3/5
syms x collect( x*(x*(x+3)+5) + 1 )
ans = x^3 + 3*x^2 + 5*x + 1
syms x y collect(x^2*y + y*x - x^2 - 2*x)
ans = (y - 1)*x^2 + (y - 2)*x
syms x y z collect ( (x+y+z) * (x-y-z) ) collect ( (x+y+z) * (x-y-z), y ) collect ( (x+y+z) * (x-y-z), z )
ans = x^2 - (y + z)^2 ans = - y^2 - 2*z*y + (x + z)*(x - z) ans = - z^2 - 2*y*z + (x + y)*(x - y)
syms x f = x^2 * exp(x); diff(f)
ans = x^2*exp(x) + 2*x*exp(x)
diff(sin(pi*x), 2)
ans = -pi^2*sin(pi*x)
syms x y g = x*y + x^3; diff(g, x) diff(g, y) diff(g, x, 2)
ans = 3*x^2 + y ans = x ans = 6*x
syms x int(-2*x/(1 + x^2)^2) syms x y z int(x*y^2 + y*z, y)
ans = 1/(x^2 + 1) ans = (y^2*(3*z + 2*x*y))/6
int(sin(x), 0 ,pi) syms x t int(2*x, [sin(t), 1])
ans = 2 ans = cos(t)^2
syms n x limit((1+x/n)^n, n, inf)
ans = exp(x)
syms x taylor(cos(x^2), x, pi, 'Order', 3)
ans = cos(pi^2) - (x - pi)^2*(sin(pi^2) + 2*pi^2*cos(pi^2)) - 2*pi*sin(pi^2)*(x - pi)
syms x y z f = x^2 + y^2 + z^2 + x*y*z; solve(f)
ans = - (y^2*z^2 - 4*y^2 - 4*z^2)^(1/2)/2 - (y*z)/2 (y^2*z^2 - 4*y^2 - 4*z^2)^(1/2)/2 - (y*z)/2
syms x y z E1 = (x - 2)^2+(y - 4)^2 - z^2; E2 = x/2 + 1 - y; [xc,yc] = solve(E1 , E2)
xc = ((4*z^2)/5 - 64/25)^(1/2) + 14/5 14/5 - ((4*z^2)/5 - 64/25)^(1/2) yc = ((4*z^2)/5 - 64/25)^(1/2)/2 + 12/5 12/5 - ((4*z^2)/5 - 64/25)^(1/2)/2