v16i065: nlmdl - Estimate nonlinear statistical models, Part03/06
Ron Gallant
arg at ccvr1.cc.ncsu.edu
Sat Jan 12 16:06:14 AEST 1991
Submitted-by: arg at ccvr1.cc.ncsu.edu (Ron Gallant)
Posting-number: Volume 16, Issue 65
Archive-name: nlmdl/part03
# This is a shell archive. Remove anything before this line, then unpack
# it by saving it into a file and typing "sh file". To overwrite existing
# files, type "sh file -c". You can also feed this as standard input via
# unshar, or by typing "sh <file", e.g.. If this archive is complete, you
# will see the following message at the end:
# "End of archive 3 (of 6)."
# Contents: ch1eg1/detail.bak ch5eg1/model.cc ch6eg2/detail.bak
# ch6eg2/hansen.dat nlmdl.opr
# Wrapped by arg at sparc on Fri Jan 11 20:14:53 1991
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f 'ch1eg1/detail.bak' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'ch1eg1/detail.bak'\"
else
echo shar: Extracting \"'ch1eg1/detail.bak'\" \(11967 characters\)
sed "s/^X//" >'ch1eg1/detail.bak' <<'END_OF_FILE'
X
X
X **********************************************************************
X * *
X * data *
X * *
X **********************************************************************
X
X
X Col 1 Col 2 Col 3 Col 4 Col 5
X
X Row 1 1.00000 0.986100 1.00000 1.00000 6.28000
X Row 2 2.00000 1.03848 0.0 1.00000 9.86000
X Row 3 3.00000 0.954820 1.00000 1.00000 9.11000
X Row 4 4.00000 1.04184 0.0 1.00000 8.43000
X Row 5 5.00000 1.02324 1.00000 1.00000 8.11000
X Row 6 6.00000 0.904750 0.0 1.00000 1.82000
X Row 7 7.00000 0.962630 1.00000 1.00000 6.58000
X Row 8 8.00000 1.05026 0.0 1.00000 5.02000
X Row 9 9.00000 0.988610 1.00000 1.00000 6.52000
X Row 10 10.0000 1.03437 0.0 1.00000 3.75000
X Row 11 11.0000 0.989820 1.00000 1.00000 9.86000
X Row 12 12.0000 1.01214 0.0 1.00000 7.31000
X Row 13 13.0000 0.667680 1.00000 1.00000 0.470000
X Row 14 14.0000 0.551070 0.0 1.00000 0.0700000
X Row 15 15.0000 0.968220 1.00000 1.00000 4.07000
X Row 16 16.0000 0.988230 0.0 1.00000 4.61000
X Row 17 17.0000 0.597590 1.00000 1.00000 0.170000
X Row 18 18.0000 0.994180 0.0 1.00000 6.99000
X Row 19 19.0000 1.01962 1.00000 1.00000 4.39000
X Row 20 20.0000 0.691630 0.0 1.00000 0.390000
X Row 21 21.0000 1.04255 1.00000 1.00000 4.73000
X Row 22 22.0000 1.04343 0.0 1.00000 9.42000
X Row 23 23.0000 0.975260 1.00000 1.00000 8.90000
X Row 24 24.0000 1.04969 0.0 1.00000 3.02000
X Row 25 25.0000 0.802190 1.00000 1.00000 0.770000
X Row 26 26.0000 1.01046 0.0 1.00000 3.31000
X Row 27 27.0000 0.951960 1.00000 1.00000 4.51000
X Row 28 28.0000 0.976580 0.0 1.00000 2.65000
X Row 29 29.0000 0.508110 1.00000 1.00000 0.0800000
X Row 30 30.0000 0.918400 0.0 1.00000 6.11000
X
X
X **********************************************************************
X * *
X * nlmdl 2.0 *
X * *
X **********************************************************************
X
X
X **********************************************************************
X * *
X * Parameter settings *
X * *
X **********************************************************************
X
X Gallant, "Nonlinear Statistical Models," Chapter 1, Figure 5, p. 35-36.
X SUR What estimation method? Code SUR, TSLS, or GMM.
X 30 Number of observations, t = 1, ..., n.
X 1 Number of equations, i.e. dimension of e.
X 0 Number of instruments, i.e. dimension of Z.
X 4 Number of parameters, i.e. dimension of theta.
X 50 Upper limit on Gauss-Newton iterations.
X 1 Number var iterates, ivar=0 means none.
X homoskedastic Code homoskedastic or heteroskedastic.
X 0 Number of moving average terms MA for var estimate.
X none Code none or Parzen, none when MA>0 is unwise.
X 1.000000e-13 Convergence tolerance, tol=1.0e-8 is reasonable.
X 1.000000e-10 Inversion tolerance, eps=1.0e-13 is reasonable
X full How much output? Code none, minimal, or full.
X
X
X **********************************************************************
X * *
X * Starting theta *
X * *
X **********************************************************************
X
X
X Col 1
X
X Row 1 -0.0486600
X Row 2 1.03884
X Row 3 -0.737920
X Row 4 -0.513620
X
Xvar_loop 0
X 1.0000000000000000e+00 1.00000000 var(1,1)
X
Xtheta_loop 0
X -4.8660000000000002e-02 -0.04866000 theta(1)
X 1.0388400000000000e+00 1.03884000 theta(2)
X -7.3792000000000002e-01 -0.73792000 theta(3)
X -5.1361999999999997e-01 -0.51362000 theta(4)
X 5.0775308202588146e-02 0.05077531 obj
X 2.4331006069174543e-02 0.02433101 D(1)
X -2.8980780513591580e-02 -0.02898078 D(2)
X -2.7779093350995165e-01 -0.27779093 D(3)
X 2.2218376750406332e-02 0.02221838 D(4)
XStep length = 1
X
Xtheta_loop 1
X -2.4328993930825459e-02 -0.02432899 theta(1)
X 1.0098592194864084e+00 1.00985922 theta(2)
X -1.0157109335099517e+00 -1.01571093 theta(3)
X -4.9140162324959363e-01 -0.49140162 theta(4)
X 3.2351520567053427e-02 0.03235152 obj
X -1.4057031119702702e-03 -0.00140570 D(1)
X 5.4557838716931487e-03 0.00545578 D(2)
X -1.0039354203311036e-01 -0.10039354 D(3)
X -1.3173231905220829e-02 -0.01317323 D(4)
XStep length = 1
X
Xtheta_loop 2
X -2.5734697042795730e-02 -0.02573470 theta(1)
X 1.0153150033581015e+00 1.01531500 theta(2)
X -1.1161044755430620e+00 -1.11610448 theta(3)
X -5.0457485515481448e-01 -0.50457486 theta(4)
X 3.0497613121941033e-02 0.03049761 obj
X -1.5509486383662779e-04 -0.00015509 D(1)
X 3.6499057445712069e-04 0.00036499 D(2)
X 4.2218191744343219e-04 0.00042218 D(3)
X -3.2672143437975117e-04 -0.00032672 D(4)
XStep length = 1
X
Xtheta_loop 3
X -2.5889791906632358e-02 -0.02588979 theta(1)
X 1.0156799939325587e+00 1.01567999 theta(2)
X -1.1156822936256185e+00 -1.11568229 theta(3)
X -5.0490157658919421e-01 -0.50490158 theta(4)
X 3.0495537203056537e-02 0.03049554 obj
X 9.7321692293949851e-08 0.00000010 D(1)
X -3.3659920289733601e-07 -0.00000034 D(2)
X -1.5374080636796911e-05 -0.00001537 D(3)
X -1.3314150623113408e-06 -0.00000133 D(4)
XStep length = 1
X
Xtheta_loop 4
X -2.5889694584940063e-02 -0.02588969 theta(1)
X 1.0156796573333557e+00 1.01567966 theta(2)
X -1.1156976677062553e+00 -1.11569767 theta(3)
X -5.0490290800425652e-01 -0.50490291 theta(4)
X 3.0495537193052542e-02 0.03049554 obj
X -3.6229294260247439e-09 -0.00000000 D(1)
X 1.2413365101829759e-08 0.00000001 D(2)
X 5.4662401885543252e-07 0.00000055 D(3)
X 4.8150321542521305e-08 0.00000005 D(4)
XStep length = 1
X
Xtheta_loop 5
X -2.5889698207869488e-02 -0.02588970 theta(1)
X 1.0156796697467207e+00 1.01567967 theta(2)
X -1.1156971210822364e+00 -1.11569712 theta(3)
X -5.0490285985393502e-01 -0.50490286 theta(4)
X 3.0495537193039934e-02 0.03049554 obj
X 1.2871068769735452e-10 0.00000000 D(1)
X -4.4115387457665704e-10 -0.00000000 D(2)
X -1.9437985890799474e-08 -0.00000002 D(3)
X -1.7124136454193738e-09 -0.00000000 D(4)
XStep length = 1
X
Xtheta_loop 6
X -2.5889698079158800e-02 -0.02588970 theta(1)
X 1.0156796693055670e+00 1.01567967 theta(2)
X -1.1156971405202223e+00 -1.11569714 theta(3)
X -5.0490286156634867e-01 -0.50490286 theta(4)
X 3.0495537193039886e-02 0.03049554 obj
X -4.5771437765728809e-12 -0.00000000 D(1)
X 1.5687617718772477e-11 0.00000000 D(2)
X 6.9121248119946195e-10 0.00000000 D(3)
X 6.0893027878551739e-11 0.00000000 D(4)
XStep length = .9
X
Xtheta_loop 7
X -2.5889698083278231e-02 -0.02588970 theta(1)
X 1.0156796693196859e+00 1.01567967 theta(2)
X -1.1156971398981310e+00 -1.11569714 theta(3)
X -5.0490286151154495e-01 -0.50490286 theta(4)
X 3.0495537193039865e-02 0.03049554 obj
X -3.1123169737974754e-13 -0.00000000 D(1)
X 1.0666118487966118e-12 0.00000000 D(2)
X 4.6999452916535867e-11 0.00000000 D(3)
X 4.1405133028926090e-12 0.00000000 D(4)
XA line search did not improve this estimate.
X
Xvar_loop 1
X 1.1729052766553795e-03 0.00117291 var(1,1)
X
Xtheta_loop 0
X -2.5889698083278231e-02 -0.02588970 theta(1)
X 1.0156796693196859e+00 1.01567967 theta(2)
X -1.1156971398981310e+00 -1.11569714 theta(3)
X -5.0490286151154495e-01 -0.50490286 theta(4)
X 2.5999999999999996e+01 26.00000000 obj
X -3.1122653130181253e-13 -0.00000000 D(1)
X 1.0666063107489048e-12 0.00000000 D(2)
X 4.6999396801498302e-11 0.00000000 D(3)
X 4.1405085967972343e-12 0.00000000 D(4)
XA line search did not improve this estimate.
X
X
X **********************************************************************
X * *
X * theta *
X * *
X **********************************************************************
X
X
X Col 1
X
X Row 1 -0.0258897
X Row 2 1.01568
X Row 3 -1.11570
X Row 4 -0.504903
X
X
X
X **********************************************************************
X * *
X * var *
X * (corrected for degrees of freedom) *
X * *
X **********************************************************************
X
X
X Col 1
X
X Row 1 0.00117291
X
X
X
X **********************************************************************
X * *
X * V *
X * (rank = 4) *
X * *
X **********************************************************************
X
X
X Col 1 Col 2 Col 3 Col 4
X
X Row 1 0.00015936 -7.872e-05 -0.00017711 -4.409e-05
X Row 2 -7.872e-05 9.876e-05 0.00060702 -1.851e-06
X Row 3 -0.00017711 0.00060702 0.0267460 0.00235621
X Row 4 -4.409e-05 -1.851e-06 0.00235621 0.00065829
X
END_OF_FILE
if test 11967 -ne `wc -c <'ch1eg1/detail.bak'`; then
echo shar: \"'ch1eg1/detail.bak'\" unpacked with wrong size!
fi
# end of 'ch1eg1/detail.bak'
fi
if test -f 'ch5eg1/model.cc' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'ch5eg1/model.cc'\"
else
echo shar: Extracting \"'ch5eg1/model.cc'\" \(8228 characters\)
sed "s/^X//" >'ch5eg1/model.cc' <<'END_OF_FILE'
X#include "model.h"
X
Xmodel::model() { }
X
Xmodel::~model() { }
X
Xvoid model::e(INTEGER t)
X{
X REAL peak, inter, base;
X REAL t1, t2, t3, t4, t5, t6, t7, t8;
X REAL r1, r2, r3;
X REAL x1, x2, x3;
X REAL y1, y2, y3;
X
X t1 = s.theta[1];
X t2 = s.theta[2];
X t3 = s.theta[3];
X t4 = s.theta[4];
X t5 = s.theta[5];
X t6 = s.theta[6];
X t7 = s.theta[7];
X t8 = s.theta[8];
X
X y1 = data.elem(1,t);
X y2 = data.elem(2,t);
X y3 = data.elem(3,t);
X
X r1 = data.elem(4,t);
X r2 = data.elem(5,t);
X r3 = data.elem(6,t);
X
X x1 = r1 - y3;
X x2 = r2 - y3;
X x3 = r3 - y3;
X
X peak = t1 + t2*x1 + t3*x2 + t4*x3;
X inter = t5 + t3*x1 + t6*x2 + t7*x3;
X base = -1 + t4*x1 + t7*x2 + t8*x3;
X
X e_t[1] = y1 - log(peak/base);
X e_t[2] = y2 - log(inter/base);
X}
X
Xvoid model::dele(INTEGER t)
X{
X REAL peak, inter, base;
X REAL t1, t2, t3, t4, t5, t6, t7, t8;
X REAL r1, r2, r3;
X REAL x1, x2, x3;
X REAL y1, y2, y3;
X
X t1 = s.theta[1];
X t2 = s.theta[2];
X t3 = s.theta[3];
X t4 = s.theta[4];
X t5 = s.theta[5];
X t6 = s.theta[6];
X t7 = s.theta[7];
X t8 = s.theta[8];
X
X y1 = data.elem(1,t);
X y2 = data.elem(2,t);
X y3 = data.elem(3,t);
X
X r1 = data.elem(4,t);
X r2 = data.elem(5,t);
X r3 = data.elem(6,t);
X
X x1 = r1 - y3;
X x2 = r2 - y3;
X x3 = r3 - y3;
X
X peak = t1 + t2*x1 + t3*x2 + t4*x3;
X inter = t5 + t3*x1 + t6*x2 + t7*x3;
X base = -1 + t4*x1 + t7*x2 + t8*x3;
X
X//e_t[1] = y1 - log(peak/base);
X//e_t[2] = y2 - log(inter/base);
X
X
X dele_t.elem(1,1) = - (1.0/peak);
X dele_t.elem(1,2) = - (1.0/peak)*x1;
X dele_t.elem(1,3) = - (1.0/peak)*x2;
X dele_t.elem(1,4) = - (1.0/peak)*x3 + (1.0/base)*x1;
X dele_t.elem(1,5) = 0;
X dele_t.elem(1,6) = 0;
X dele_t.elem(1,7) = (1.0/base)*x2;
X dele_t.elem(1,8) = (1.0/base)*x3;
X
X dele_t.elem(2,1) = 0;
X dele_t.elem(2,2) = 0;
X dele_t.elem(2,3) = - (1.0/inter)*x1;
X dele_t.elem(2,4) = (1.0/base)*x1;
X dele_t.elem(2,5) = - (1.0/inter);
X dele_t.elem(2,6) = - (1.0/inter)*x2;
X dele_t.elem(2,7) = - (1.0/inter)*x3 + (1.0/base)*x2;
X dele_t.elem(2,8) = (1.0/base)*x3;
X}
X
Xint model::initialize()
X{
X if (s.p != 8 || s.M != 2) {cerr << "Error, model::initialize.\n"; exit(1);}
X
X#ifdef GNU_GPP_COMPILER
X
X#ifdef USE_ATT_STYLE_IO_WITH_GNU
X
X filebuf electa_buf;
X if( electa_buf.open("electa.dat", input) == 0 ){
X cerr << "Error, model::model, Cannot open electa.dat\n";
X exit(1);
X }
X istream electa(&electa_buf);
X
X
X filebuf electc1_buf;
X if( electc1_buf.open("electc1.dat", input) == 0 ){
X cerr << "Error, model::model, Cannot open electc1.dat\n";
X exit(1);
X }
X istream electc1(&electc1_buf);
X
X filebuf electc2_buf;
X if( electc2_buf.open("electc2.dat", input) == 0 ){
X cerr << "Error, model::model, Cannot open electc2.dat\n";
X exit(1);
X }
X istream electc2(&electc2_buf);
X
X#endif
X
X#ifdef USE_GNU_STYLE_IO_WITH_GNU
X
X istream electa("electa.dat",io_readonly,a_useonly);
X if (!electa) {
X cerr << "Error, model::model, Cannot open electa.dat\n";
X exit(1);
X }
X
X istream electc1("electc1.dat",io_readonly,a_useonly);
X if (!electc1) {
X cerr << "Error, model::model, Cannot open electc1.dat\n";
X exit(1);
X }
X
X istream electc2("electc2.dat",io_readonly,a_useonly);
X if (!electc2) {
X cerr << "Error, model::model, Cannot open electc2.dat\n";
X exit(1);
X }
X
X#endif
X
X#endif
X
X#ifdef TURBO_CPP_COMPILER
X
X ifstream electa("electa.dat");
X if (!electa) {
X cerr << "Error, model::model, Cannot open electa.dat\n";
X exit(1);
X }
X
X ifstream electc1("electc1.dat");
X if (!electc1) {
X cerr << "Error, model::model, Cannot open electc1.dat\n";
X exit(1);
X }
X
X ifstream electc2("electc2.dat");
X if (!electc2) {
X cerr << "Error, model::model, Cannot open electc2.dat\n";
X exit(1);
X }
X#endif
X
X double t1, treat, base, inter, peak, expend;
X double t2, famsize, income, sqfeet, heatlos, range, wash, dry, cac, wac;
X double t3, single, duplex, mobile, hwh, frez, ref;
X
X data.resize(20,224);
X
X int missing;
X
X s.n = 224;
X
X INTEGER t0 = 1;
X
X for (INTEGER t=1; t<=224; t++) {
X
X missing = 0;
X
X electa >> t1;
X electa >> treat;
X electa >> base;
X electa >> inter;
X electa >> peak;
X electa >> expend;
X
X electc1 >> t2;
X electc1 >> famsize; if(famsize == -9999.) missing = 1; //used for d8 d10
X electc1 >> income; if(income == -9999.) missing = 1; //used for d0
X electc1 >> sqfeet; if(sqfeet == -9999.) missing = 1; //used for d3
X electc1 >> heatlos; if(heatlos == -9999.) missing = 1; //used for d6
X electc1 >> range; if(range == -9999.) missing = 1; //used for d13
X electc1 >> wash; if(wash == -9999.) missing = 1; //used for d9
X electc1 >> dry; if(dry == -9999.) missing = 1; //used for d10
X electc1 >> cac; if(cac == -9999.) missing = 1; //used for d6
X electc1 >> wac; if(wac == -9999.) missing = 1; //used for d7
X
X electc2 >> t3;
X electc2 >> single; // if(single == -9999.) missing = 1; //not used
X electc2 >> duplex; if(duplex == -9999.) missing = 1; //used for d4
X electc2 >> mobile; if(mobile == -9999.) missing = 1; //used for d5
X electc2 >> hwh; if(hwh == -9999.) missing = 1; //used for d8 d9
X electc2 >> frez; if(frez == -9999.) missing = 1; //used for d12
X electc2 >> ref; if(ref == -9999.) missing = 1; //used for d11
X
X if (t1 != t2 || t1 != t3) {
X cerr << "Error, model::initialize, t, t1, t2, t3 = "
X << t << " " << t1 << " " << t2 << " " << t3 << ".\n";
X exit(1);
X }
X
X double p1, p2, p3;
X
X p1 = p2 = p3 = 0;
X
X if (treat == 1) { p1=3.90; p2=2.86; p3=1.06; }
X if (treat == 2) { p1=3.90; p2=2.86; p3=1.78; }
X if (treat == 3) { p1=3.90; p2=3.90; p3=1.06; }
X if (treat == 4) { p1=3.90; p2=3.90; p3=1.78; }
X if (treat == 5) { p1=5.06; p2=3.34; p3=1.37; }
X if (treat == 6) { p1=6.56; p2=2.86; p3=1.06; }
X if (treat == 7) { p1=6.56; p2=2.86; p3=1.78; }
X if (treat == 8) { p1=6.56; p2=3.90; p3=1.06; }
X if (treat == 9) { p1=6.56; p2=3.90; p3=1.78; }
X
X
X double y1 = log(peak/base);
X double y2 = log(inter/base);
X double y3 = log(expend);
X
X double r1 = log(p1);
X double r2 = log(p2);
X double r3 = log(p3);
X
X double d0, d1, d2, d3, d4, d5, d6, d7, d8, d9, d10, d11, d12, d13;
X
X if (missing) {
X cout << "Missing observation at t = " << t << ".\n";
X d0= d1= d2= d3= d4= d5= d6= d7= d8= d9= d10= d11= d12= d13= -9999.;
X }
X else {
X d0 = 1;
X d1 = log(10*p1+6*p2+8*p3);
X d2 = log(income);
X d3 = log(sqfeet);
X d4 = duplex;
X d5 = mobile;
X d6 = cac*log(heatlos);
X d7 = 0; if (wac>0) d7=log(wac);
X d8 = 0; if (hwh>0) d8=log(famsize+1);
X d9 = 0; if (hwh>0 && wash>0) d9=1;
X d10 = 0; if (dry>0) d10=log(famsize+1);
X d11 = 0; if (ref>0) d11=log(ref);
X d12 = 0; if (frez>0) d12=log(frez);
X d13 = range;
X }
X
X if (strcmp(s.method,"SUR") == 0 || !missing ) {
X data.elem(1,t0) = (REAL)y1;
X data.elem(2,t0) = (REAL)y2;
X data.elem(3,t0) = (REAL)y3;
X data.elem(4,t0) = (REAL)r1;
X data.elem(5,t0) = (REAL)r2;
X data.elem(6,t0) = (REAL)r3;
X data.elem(7,t0) = (REAL)d0;
X data.elem(8,t0) = (REAL)d1;
X data.elem(9,t0) = (REAL)d2;
X data.elem(10,t0) = (REAL)d3;
X data.elem(11,t0) = (REAL)d4;
X data.elem(12,t0) = (REAL)d5;
X data.elem(13,t0) = (REAL)d6;
X data.elem(14,t0) = (REAL)d7;
X data.elem(15,t0) = (REAL)d8;
X data.elem(16,t0) = (REAL)d9;
X data.elem(17,t0) = (REAL)d10;
X data.elem(18,t0) = (REAL)d11;
X data.elem(19,t0) = (REAL)d12;
X data.elem(20,t0) = (REAL)d13;
X t0++;
X }
X else {
X s.n--;
X }
X }
X
X s.starting="tmp.dat";
X s.to(s.starting);
X
X// if (strcmp(s.detail,"full") == 0) cout << starbox("/data//_") << data;
X
X return 1;
X}
X
Xvoid model::Z(INTEGER t)
X{
X INTEGER i;
X
X if (s.K < 4 || s.K > 17) {cerr << "Error, model::Z.\n"; exit(1);}
X
X for (i = 1; i<=3; i++) Z_t[i] = data.elem(i+3,t);
X for (i = 4; i<=s.K; i++) Z_t[i] = data.elem(i+3,t);
X
X}
X
Xint model::terminate() { return 0; }
END_OF_FILE
if test 8228 -ne `wc -c <'ch5eg1/model.cc'`; then
echo shar: \"'ch5eg1/model.cc'\" unpacked with wrong size!
fi
# end of 'ch5eg1/model.cc'
fi
if test -f 'ch6eg2/detail.bak' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'ch6eg2/detail.bak'\"
else
echo shar: Extracting \"'ch6eg2/detail.bak'\" \(7346 characters\)
sed "s/^X//" >'ch6eg2/detail.bak' <<'END_OF_FILE'
X
X
X **********************************************************************
X * *
X * nlmdl 2.0 *
X * *
X **********************************************************************
X
X
X **********************************************************************
X * *
X * Parameter settings *
X * *
X **********************************************************************
X
X Gallant, "Nonlinear Statistical Models," Chapter 6, Figure 2, p. 449-450.
X GMM What estimation method? Code SUR, TSLS, or GMM.
X 238 Number of observations, t = 1, ..., n.
X 1 Number of equations, i.e. dimension of e.
X 3 Number of instruments, i.e. dimension of Z.
X 2 Number of parameters, i.e. dimension of theta.
X 20 Upper limit on Gauss-Newton iterations.
X 1 Number var iterates, ivar=0 means none.
X homoskedastic Code homoskedastic or heteroskedastic.
X 0 Number of moving average terms MA for var estimate.
X Parzen Code none or Parzen, none when MA>0 is unwise.
X 1.000000e-05 Convergence tolerance, tol=1.0e-8 is reasonable.
X 1.000000e-13 Inversion tolerance, eps=1.0e-13 is reasonable
X full How much output? Code none, minimal, or full.
X
X
X **********************************************************************
X * *
X * Starting theta *
X * *
X **********************************************************************
X
X
X Col 1
X
X Row 1 -0.400000
X Row 2 0.900000
X
Xvar_loop 0
X 2.3800000000000000e+02 238.00000000 var(1,1)
X 2.3848513503508150e+02 238.48513504 var(2,1)
X 2.3866807009536583e+02 238.66807010 var(3,1)
X 2.3848513503508150e+02 238.48513504 var(1,2)
X 2.3897620940526255e+02 238.97620941 var(2,2)
X 2.3916481122626234e+02 239.16481123 var(3,2)
X 2.3866807009536583e+02 238.66807010 var(1,3)
X 2.3916481122626234e+02 239.16481123 var(2,3)
X 2.3975964092305995e+02 239.75964092 var(3,3)
X
Xtheta_loop 0
X -4.0000000000000002e-01 -0.40000000 theta(1)
X 9.0000000000000002e-01 0.90000000 theta(2)
X 2.2977593855825589e+00 2.29775939 obj
X -4.9910043429446205e-01 -0.49910043 D(1)
X 9.8933364963875864e-02 0.09893336 D(2)
XStep length = 1
X
Xtheta_loop 1
X -8.9910043429446207e-01 -0.89910043 theta(1)
X 9.9893336496387586e-01 0.99893336 theta(2)
X 2.7751760221041621e-03 0.00277518 obj
X 5.0366043831498272e-02 0.05036604 D(1)
X -4.4414610501568677e-06 -0.00000444 D(2)
XStep length = 1
X
Xtheta_loop 2
X -8.4873439046296384e-01 -0.84873439 theta(1)
X 9.9892892350282569e-01 0.99892892 theta(2)
X 2.7717868638761308e-03 0.00277179 obj
X -1.1734797022229638e-04 -0.00011735 D(1)
X 2.1096605255709332e-07 0.00000021 D(2)
XStep length = 1
X
Xtheta_loop 3
X -8.4885173843318618e-01 -0.84885174 theta(1)
X 9.9892913446887821e-01 0.99892913 theta(2)
X 2.7717868580610532e-03 0.00277179 obj
X 2.6798358348564710e-07 0.00000027 D(1)
X -5.5350439549602871e-10 -0.00000000 D(2)
XTolerence check passed.
X
Xvar_loop 1
X 4.0582168346665826e-01 0.40582168 var(1,1)
X 4.0643410965960231e-01 0.40643411 var(2,1)
X 3.9873658267692480e-01 0.39873658 var(3,1)
X 4.0643410965960231e-01 0.40643411 var(1,2)
X 4.0705487414860658e-01 0.40705487 var(2,2)
X 3.9936326231932223e-01 0.39936326 var(3,2)
X 3.9873658267692480e-01 0.39873658 var(1,3)
X 3.9936326231932223e-01 0.39936326 var(2,3)
X 3.9272327136336027e-01 0.39272327 var(3,3)
X
Xtheta_loop 0
X -8.4885173843318618e-01 -0.84885174 theta(1)
X 9.9892913446887821e-01 0.99892913 theta(2)
X 1.2467135787066870e+00 1.24671358 obj
X -1.8474288819635265e-01 -0.18474289 D(1)
X -6.7258036762946278e-04 -0.00067258 D(2)
XStep length = 1
X
Xtheta_loop 1
X -1.0335946266295388e+00 -1.03359463 theta(1)
X 9.9825655410124869e-01 0.99825655 theta(2)
X 1.0569217932960528e+00 1.05692179 obj
X 7.2887425999781236e-05 0.00007289 D(1)
X -8.3045884673608490e-07 -0.00000083 D(2)
XStep length = 1
X
Xtheta_loop 2
X -1.0335217392035390e+00 -1.03352174 theta(1)
X 9.9825572364240200e-01 0.99825572 theta(2)
X 1.0569217148505041e+00 1.05692171 obj
X -9.1931859179544446e-08 -0.00000009 D(1)
X 1.8525999753608607e-10 0.00000000 D(2)
XTolerence check passed.
X
X
X **********************************************************************
X * *
X * theta *
X * *
X **********************************************************************
X
X
X Col 1
X
X Row 1 -1.03352
X Row 2 0.998256
X
X
X
X **********************************************************************
X * *
X * var *
X * (no degrees of freedom corrections) *
X * *
X **********************************************************************
X
X
X Col 1 Col 2 Col 3
X
X Row 1 0.405822 0.406434 0.398737
X Row 2 0.406434 0.407055 0.399363
X Row 3 0.398737 0.399363 0.392723
X
X
X
X **********************************************************************
X * *
X * V *
X * (rank = 2) *
X * *
X **********************************************************************
X
X
X Col 1 Col 2
X
X Row 1 3.58009 -0.00721268
X Row 2 -0.00721268 2.060e-05
X
END_OF_FILE
if test 7346 -ne `wc -c <'ch6eg2/detail.bak'`; then
echo shar: \"'ch6eg2/detail.bak'\" unpacked with wrong size!
fi
# end of 'ch6eg2/detail.bak'
fi
if test -f 'ch6eg2/hansen.dat' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'ch6eg2/hansen.dat'\"
else
echo shar: Extracting \"'ch6eg2/hansen.dat'\" \(11280 characters\)
sed "s/^X//" >'ch6eg2/hansen.dat' <<'END_OF_FILE'
X 381.9 176.6850 0.0093695102 0.6818539
X 383.7 176.9050 0.0093310997 0.6823039
X 388.3 177.1460 0.0049904501 0.6814319
X 385.5 177.3650 0.0383739690 0.6830091
X 389.7 177.5910 0.0204769890 0.6846292
X 390.0 177.8300 0.0007165600 0.6876923
X 389.2 178.1010 0.0371922290 0.6893628
X 390.7 178.3760 -0.0113433900 0.6910673
X 393.6 178.6570 -0.0472779090 0.6930894
X 394.2 178.9210 0.0164727200 0.6945713
X 394.1 179.1530 0.0194594210 0.6950013
X 396.5 179.3860 0.0296911900 0.6958386
X 396.8 179.5970 -0.0664901060 0.6960685
X 395.4 179.7880 0.0114439700 0.6967628
X 399.1 180.0070 -0.0114419700 0.6983212
X 404.2 180.2220 -0.0163223000 0.7013855
X 399.8 180.4440 0.0328373610 0.7016008
X 401.3 180.6710 0.0231378990 0.7024670
X 402.0 180.9450 -0.0210754290 0.7034826
X 400.4 181.2380 0.0296860300 0.7047952
X 400.2 181.5280 -0.0568203400 0.7061469
X 402.9 181.7960 -0.0045937700 0.7078680
X 403.8 182.0420 0.0472565590 0.7100049
X 401.6 182.2870 0.0478186380 0.7109064
X 404.0 182.5200 0.0654135870 0.7106436
X 405.7 182.7420 0.0364446900 0.7108701
X 409.4 182.9920 0.0318523910 0.7105520
X 410.1 183.2170 0.0058811302 0.7100707
X 412.1 183.4520 0.0251120610 0.7100218
X 412.4 183.6910 -0.0296279510 0.7109602
X 410.4 183.9580 0.0303546690 0.7127193
X 411.5 184.2430 0.0251584590 0.7132442
X 413.7 184.5240 -0.0189705600 0.7138023
X 415.9 184.7830 0.0265103900 0.7136331
X 419.0 185.0160 0.0470347110 0.7136038
X 420.5 185.2420 0.0006585300 0.7143876
X 420.8 185.4520 -0.0358958100 0.7155418
X 420.6 185.6500 0.0197925490 0.7177841
X 423.4 185.8740 -0.0055647301 0.7191781
X 424.8 186.0870 -0.0615112410 0.7203390
X 427.0 186.3140 -0.0834698080 0.7206089
X 425.2 186.5380 -0.0809235570 0.7208373
X 427.0 186.7900 0.0659098630 0.7203747
X 428.5 187.0580 0.0226466300 0.7218203
X 431.8 187.3230 -0.0487761200 0.7264938
X 431.0 187.5740 0.0039394898 0.7262181
X 433.6 187.7960 0.1114552000 0.7269373
X 434.1 188.0130 0.0139081300 0.7270214
X 434.7 188.2130 0.0508059190 0.7292386
X 433.7 188.3870 -0.0226500000 0.7299977
X 436.2 188.5800 0.0347222910 0.7297111
X 437.0 188.7900 0.0479895880 0.7295194
X 436.9 189.0180 0.0206833590 0.7308308
X 440.2 189.2420 -0.0178168900 0.7319400
X 442.1 189.4960 -0.0018435300 0.7335444
X 445.6 189.7610 0.0536292610 0.7349641
X 443.8 190.0280 -0.0126571200 0.7347905
X 444.2 190.2650 0.0286585090 0.7361549
X 445.8 190.4720 -0.0047020698 0.7375505
X 449.5 190.6680 0.0221940800 0.7385984
X 450.1 190.8580 0.0256042290 0.7398356
X 453.7 191.0470 0.0181333610 0.7399162
X 456.6 191.2450 0.0173465290 0.7402540
X 456.7 191.4470 0.0051271599 0.7407488
X 462.1 191.6660 0.0166149310 0.7405324
X 463.8 191.8890 0.0158939310 0.7416990
X 466.0 192.1310 0.0202899800 0.7429185
X 468.5 192.3760 -0.0109651800 0.7432231
X 468.0 192.6310 0.0315313120 0.7448718
X 470.0 192.8470 0.0100951200 0.7457447
X 468.0 193.0390 0.0013465700 0.7465812
X 474.4 193.2230 0.0034312201 0.7476813
X 474.5 193.3930 0.0377800500 0.7481560
X 477.4 193.5400 0.0066147698 0.7488479
X 474.5 193.7090 -0.0107356600 0.7506849
X 479.6 193.8880 0.0345496910 0.7525021
X 481.2 194.0870 -0.0047443998 0.7554032
X 479.5 194.3030 -0.0505878400 0.7593326
X 484.3 194.5280 0.0169978100 0.7602726
X 485.3 194.7610 0.0299301090 0.7601484
X 488.7 194.9970 0.0323472920 0.7605893
X 497.2 195.1950 0.0293272190 0.7626710
X 497.1 195.3720 0.0008636100 0.7648361
X 499.0 195.5390 0.0121703600 0.7671343
X 500.1 195.6880 0.0100357400 0.7696461
X 501.5 195.8310 -0.0102875900 0.7730808
X 502.9 195.9990 -0.0215729900 0.7757009
X 505.8 196.1780 0.0233628400 0.7785686
X 504.8 196.3720 -0.0509349700 0.7793185
X 507.5 196.5600 -0.0109703900 0.7812808
X 510.9 196.7620 -0.0118703500 0.7827363
X 508.3 196.9840 -0.0748946070 0.7867401
X 510.2 197.2070 -0.0066132201 0.7894943
X 509.8 197.3980 0.0464050400 0.7910945
X 512.1 197.5720 0.0138342800 0.7922281
X 513.5 197.7360 0.0047225100 0.7933788
X 516.0 197.8920 0.0838221310 0.7941861
X 517.7 198.0370 0.0098125497 0.7948619
X 519.0 198.2060 0.0433843100 0.7959538
X 521.1 198.3630 0.0420965220 0.7965842
X 521.0 198.5370 -0.0415207000 0.7988484
X 523.1 198.7120 0.0232013710 0.8015676
X 522.1 198.9110 0.0482556600 0.8038690
X 525.5 199.1130 -0.0056581302 0.8058991
X 528.2 199.3110 0.0336121990 0.8076486
X 524.9 199.4980 -0.0276739710 0.8094875
X 527.9 199.6570 0.0078005102 0.8124645
X 531.9 199.8080 0.0307225010 0.8155668
X 533.0 199.9200 -0.0389530290 0.8200750
X 533.9 200.0560 -0.0311505910 0.8231879
X 539.8 200.2080 0.0068653398 0.8262319
X 540.0 200.3610 0.0898963210 0.8285185
X 541.2 200.5360 0.0230161700 0.8318551
X 547.8 200.7060 0.0118528200 0.8335159
X 550.9 200.8980 -0.0211507900 0.8359049
X 552.4 201.0950 0.0163246690 0.8388849
X 551.0 201.2900 0.0422958400 0.8421053
X 552.1 201.4660 0.0111537600 0.8462235
X 556.7 201.6210 0.0562853110 0.8492905
X 554.1 201.7600 -0.0372401590 0.8521928
X 557.0 201.8810 -0.0072337599 0.8560144
X 561.2 202.0230 -0.0502119700 0.8578047
X 560.6 202.1610 0.0314719300 0.8619336
X 561.9 202.3310 0.0213753300 0.8667023
X 566.5 202.5070 0.0029275999 0.8704325
X 563.9 202.6770 -0.0623450020 0.8751552
X 565.9 202.8770 -0.0630705430 0.8784238
X 569.4 203.0900 0.0504970810 0.8816298
X 568.2 203.3020 -0.0220447110 0.8856037
X 573.1 203.5000 0.0547974710 0.8888501
X 572.5 203.6750 -0.0314589110 0.8939738
X 572.4 203.8490 -0.0180749090 0.8981481
X 577.2 204.0080 -0.0763489600 0.9019404
X 578.1 204.1560 0.0597185420 0.9058986
X 577.7 204.3350 -0.0023485899 0.9077376
X 577.1 204.5050 -0.0995834470 0.9123202
X 580.3 204.6920 -0.0611347710 0.9155609
X 582.0 204.8780 -0.0502832790 0.9176976
X 582.8 205.0860 0.0746088620 0.9208991
X 584.7 205.2940 0.0502020900 0.9235505
X 588.5 205.5070 0.0426676610 0.9276126
X 587.3 205.7070 -0.0160981810 0.9324025
X 587.6 205.8840 0.0521828200 0.9361811
X 592.6 206.0760 0.0617985200 0.9399258
X 592.2 206.2420 0.0492740680 0.9414049
X 594.5 206.3930 0.0149685600 0.9434819
X 592.4 206.5670 0.0441647210 0.9469953
X 596.1 206.7260 0.0341992900 0.9506794
X 596.3 206.8910 -0.0365711710 0.9548885
X 598.5 207.0530 0.0043891501 0.9597327
X 597.3 207.2370 -0.0398038400 0.9630002
X 599.1 207.4330 0.0409017500 0.9679519
X 601.1 207.6270 -0.0056930701 0.9698885
X 601.7 207.8000 -0.0395274310 0.9729101
X 604.9 207.9490 -0.0000956400 0.9752025
X 608.8 208.0880 0.0907427070 0.9797963
X 607.9 208.1960 0.0241155400 0.9825629
X 610.3 208.3100 0.0308808110 0.9875471
X 618.9 208.4470 0.0091922097 0.9891743
X 620.6 208.5690 0.0066767102 0.9911376
X 622.3 208.7120 0.0176741590 0.9942150
X 623.7 208.8460 -0.0221355410 0.9961520
X 627.6 208.9880 -0.0018799200 0.9996813
X 629.7 209.1530 0.0382015890 1.0031760
X 631.7 209.3170 -0.0064313002 1.0079150
X 638.2 209.4570 0.0099495398 1.0117520
X 639.8 209.5840 0.0497901590 1.0151610
X 640.7 209.7110 0.0116306100 1.0190420
X 643.4 209.8090 -0.0253177100 1.0247120
X 645.3 209.9050 -0.0398146990 1.0309930
X 643.3 210.0340 -0.0054550399 1.0399500
X 642.1 210.1540 -0.0464594700 1.0478120
X 643.2 210.2860 -0.0183557910 1.0541040
X 646.0 210.4100 -0.0088413004 1.0603720
X 651.9 210.5560 0.0521348010 1.0632000
X 643.4 210.7150 -0.0302029100 1.0792660
X 651.3 210.8630 0.0522540810 1.0815290
X 649.5 210.9840 -0.0018884100 1.0896070
X 651.3 211.0970 -0.1165516000 1.0993400
X 647.7 211.2070 0.0153318600 1.1093100
X 648.4 211.3110 -0.0013036400 1.1215300
X 646.2 211.4110 0.0038444500 1.1363350
X 645.9 211.5220 -0.0243075400 1.1489390
X 648.6 211.6370 -0.0433935780 1.1558740
X 649.3 211.7720 -0.0352610800 1.1667950
X 650.3 211.9010 -0.0193944290 1.1737660
X 653.5 212.0510 -0.0730255170 1.1802600
X 654.5 212.2160 -0.0852853730 1.1926660
X 652.7 212.3830 -0.1098341000 1.2043820
X 654.5 212.5180 0.1671594000 1.2122230
X 651.2 212.6370 -0.0397416390 1.2205160
X 650.3 212.7480 -0.0234328400 1.2278950
X 653.7 212.8440 0.1358016000 1.2337460
X 657.4 212.9390 0.0607054380 1.2376030
X 659.4 213.0560 0.0293416310 1.2406730
X 659.7 213.1870 0.0470072290 1.2457180
X 670.4 213.3930 0.0546782990 1.2502980
X 669.7 213.5590 0.0517648310 1.2593700
X 668.3 213.7410 -0.0637501480 1.2721830
X 670.1 213.9000 -0.0203062710 1.2786150
X 670.2 214.0550 -0.0366309580 1.2821550
X 670.8 214.2000 0.0609995690 1.2904000
X 674.1 214.3210 0.0314961600 1.2966920
X 677.4 214.4460 -0.0105694800 1.3039560
X 684.3 214.5610 0.1251743000 1.3081980
X 682.9 214.6550 0.0012425600 1.3069260
X 687.1 214.7620 0.0300192200 1.3092710
X 690.6 214.8810 -0.0108725300 1.3132060
X 688.7 215.0180 -0.0088088503 1.3206040
X 695.0 215.1520 0.0472505990 1.3256120
X 696.8 215.3110 -0.0073758499 1.3307980
X 699.6 215.4780 0.0005799900 1.3381930
X 702.5 215.6420 0.0261333100 1.3449110
X 705.6 215.7920 -0.0214380810 1.3520410
X 709.7 215.9240 0.0046152598 1.3587430
X 715.8 216.0570 0.0585772800 1.3657450
X 717.6 216.1860 -0.0398427810 1.3720740
X 719.3 216.3000 -0.0162227190 1.3831500
X 716.5 216.4360 -0.0106509200 1.3884160
X 719.1 216.5650 0.0038957901 1.3953550
X 722.6 216.7120 -0.0126387400 1.4014670
X 721.5 216.8630 0.0509454310 1.4105340
X 728.3 217.0300 -0.0156951400 1.4159000
X 727.0 217.2070 -0.0140849800 1.4244840
X 729.1 217.3740 0.0006794800 1.4295710
X 735.7 217.5230 -0.0394544790 1.4350960
X 739.4 217.6550 0.0419719890 1.4442790
X 740.1 217.7850 0.0052549900 1.4508850
X 738.0 217.8810 -0.0568409600 1.4581300
X 744.8 217.9870 -0.0121089800 1.4663000
X 750.5 218.1310 0.0318689010 1.4743500
X 750.4 218.2610 0.0833722430 1.4862740
X 750.3 218.4040 0.0186665390 1.5033990
X 753.1 218.5480 -0.0129163500 1.5146730
X 755.6 218.7200 0.0564879100 1.5199840
X 761.1 218.9090 0.0372171590 1.5284460
X 765.4 219.0780 -0.0063229799 1.5412860
X 765.2 219.2360 -0.1017461000 1.5541030
X 768.0 219.3840 0.0313147900 1.5640620
X 774.1 219.5300 0.0166718100 1.5694350
END_OF_FILE
if test 11280 -ne `wc -c <'ch6eg2/hansen.dat'`; then
echo shar: \"'ch6eg2/hansen.dat'\" unpacked with wrong size!
fi
# end of 'ch6eg2/hansen.dat'
fi
if test -f 'nlmdl.opr' -a "${1}" != "-c" ; then
echo shar: Will not clobber existing file \"'nlmdl.opr'\"
else
echo shar: Extracting \"'nlmdl.opr'\" \(9187 characters\)
sed "s/^X//" >'nlmdl.opr' <<'END_OF_FILE'
X/* ---------------------------------------------------------------------------
X
Xnlmdl: nlndl.opr
X
Xnlmdl is a C++ implementation of the statistical methods in A. Ronald
XGallant, "Nonlinear Statistical Models," New York: John Wiley and Sons,
X1987, ISBN 0-471-80260-3, using a matrix class realmat that is distributed
Xwith it. The header files nlmdl.h and realmat.h describe the use of the
Xprogram and matrix class, respectively.
X
XCopyright (C) 1990.
X
XA. Ronald Gallant
XP.O. Box 5513
XRaleigh NC 27650-5513
XUSA
X
XPermission to use, copy, modify, and distribute this software and its
Xdocumentation for any purpose and without fee is hereby granted, provided
Xthat the above copyright notice appear in all copies and that both that
Xcopyright notice and this permission notice appear in supporting
Xdocumentation.
X
XThis software is provided "as is" without any expressed or implied warranty.
X
X------------------------------------------------------------------------------
X
XThis is the collection of operators used by nlmdl. They act on s, of class
Xstatus, using m, of class model; s must have been filled in and the function
Xinitialize(), of class model, must have been called before using any of these
Xoperators.
X
X--------------------------------------------------------------------------- */
X
Xint line_search(char* msg)
X{
X REAL norm_D, norm_theta;
X REAL step_length;
X REAL save_obj = s.obj;
X realmat save_theta = s.theta;
X realmat ss;
X
X ss = T(s.D) * s.D;
X norm_D = sqrt(ss[1]);
X
X ss = T(s.theta) * s.theta;
X norm_theta = sqrt(ss[1]);
X
X if ( norm_D < s.tol*(norm_theta + s.eps) ) {
X strcpy(msg,"Tolerence check passed.\n");
X return 1;
X }
X
X for (step_length=1; step_length >= 0.5; step_length -= 0.1) {
X
X s.theta = save_theta + step_length * s.D;
X
X (*opr_obj)();
X
X if (s.obj < save_obj) {
X sprintf(msg,"Step length = %g \n",step_length);
X return 0 ;
X }
X }
X
X for (step_length=0.25; save_theta != s.theta; step_length *= (REAL)0.25) {
X
X s.theta = save_theta + step_length * s.D;
X
X (*opr_obj)();
X
X if (s.obj < save_obj) {
X sprintf(msg,"Step length = %g \n",step_length);
X return 0;
X }
X }
X
X strcpy(msg,"A line search did not improve this estimate.\n");
X s.obj = save_obj;
X s.theta = save_theta;
X return 1;
X}
X
X
Xvoid SUR_obj()
X{
X INTEGER t;
X realmat sse(1,1,REAL_ZERO );
X realmat varinv = invpsd(s.var,s.eps);
X
X for (t=1; t<=s.n; t++) {
X m.e(t);
X sse += T(m.e_t) * varinv * m.e_t;
X }
X s.obj = sse[1];
X}
X
X
Xvoid SUR_mgn()
X{
X INTEGER t;
X realmat sse(1,1,REAL_ZERO );
X realmat varinv = invpsd(s.var,s.eps);
X
X s.V.resize(s.p, s.p, REAL_ZERO );
X s.D.resize(s.p, 1, REAL_ZERO );
X
X for(t=1; t<=s.n; t++) {
X m.e(t);
X m.dele(t);
X sse += T(m.e_t) * varinv * m.e_t;
X s.V += T(m.dele_t) * varinv * m.dele_t;
X s.D += T(m.dele_t) * varinv * m.e_t;
X }
X
X s.obj = sse[1];
X s.V = invpsd(s.V,s.eps);
X s.D = - s.V * s.D;
X}
X
X
Xvoid SUR_var(int var_loop)
X{
X INTEGER i,j,t;
X REAL fn;
X
X if (rows(s.var)==0 || cols(s.var)==0) {
X s.var.resize(s.M,s.M,REAL_ZERO );
X for (i=1; i<=s.M; i++) s.var.elem(i,i)=1.0;
X }
X
X if (var_loop == 0) return;
X
X for (i=1; i<=s.M; i++)
X for (j=1; j<=s.M; j++)
X s.var.elem(i,j) = REAL_ZERO ;
X
X for (t=1; t<=s.n; t++) {
X m.e(t);
X s.var += m.e_t * T(m.e_t);
X }
X
X if (s.M == 1 && s.n > s.p) {
X fn = s.n - s.p;
X strcpy(s.df,"corrected");
X }
X else {
X fn = s.n;
X strcpy(s.df,"uncorrected");
X }
X
X for (i=1; i<=s.M; i++)
X for (j=1; j<=s.M; j++)
X s.var.elem(i,j) = s.var.elem(i,j)/fn;
X}
X
X
Xvoid SUR_V()
X{
X INTEGER t;
X INTEGER tau;
X INTEGER i;
X REAL x,weight;
X realmat varinv = invpsd(s.var,s.eps);
X
X if (strcmp(s.vartype,"heteroskedastic") == 0 || s.MA > 0) {
X
X realmat I(s.p,s.p,REAL_ZERO );
X realmat I_tau;
X realmat e_t_lag, dele_t_lag;
X
X
X for (tau=0; tau<=s.MA; tau++) {
X
X I_tau.resize(s.p,s.p,REAL_ZERO );
X
X for (t=tau+1; t<=s.n; t++) {
X m.e(t-tau);
X e_t_lag = m.e_t;
X
X m.dele(t-tau);
X dele_t_lag = m.dele_t;
X
X m.e(t);
X m.dele(t);
X
X I_tau +=
X T(m.dele_t) * varinv * m.e_t * T(e_t_lag) * varinv * dele_t_lag;
X }
X
X
X if (strcmp(s.weights,"Parzen") == 0 && tau > 0) {
X x = tau/(REAL)s.MA;
X
X if ( x < 0.5 )
X weight = 1.0 - 6.0*pow(x,2) + 6.0*pow(x,3);
X else
X weight = 2.0*pow((1.0 - x),3);
X }
X else {
X weight = 1.0;
X }
X
X I = I + weight*I_tau;
X
X if (tau > 0)
X I = I + weight*T(I_tau);
X }
X
X s.V = s.V * I * s.V;
X }
X
X s.rank = 0;
X for (i=1; i<=s.p; i++)
X if (s.V.elem(i,i) > REAL_ZERO ) s.rank++;
X
X return;
X}
X
X
Xrealmat qZ;
X
Xvoid opr_qZ(INTEGER t)
X{
X INTEGER alpha, i, ii;
X
X m.e(t);
X m.Z(t);
X
X for (alpha=1; alpha<=s.M; alpha++) {
X for (i=1; i<=s.K; i++) {
X
X ii = s.K*(alpha-1) + i;
X
X qZ[ii] = m.e_t[alpha] * m.Z_t[i];
X }
X }
X}
X
X
Xrealmat QZ;
X
Xvoid opr_QZ(INTEGER t)
X{
X INTEGER alpha, i, ii, k;
X
X m.dele(t);
X m.Z(t);
X
X for (k=1; k<=s.p; k++) {
X for (alpha=1; alpha<=s.M; alpha++) {
X for (i=1; i<=s.K; i++) {
X
X ii = s.K*(alpha-1) + i;
X
X QZ.elem(ii,k) = m.dele_t.elem(alpha,k) * m.Z_t[i];
X }
X }
X }
X}
X
X
Xrealmat qZ_mts;
X
Xvoid opr_qZ_mts()
X{
X qZ_mts.resize(s.M*s.K, 1, REAL_ZERO );
X
X for (INTEGER t=1; t<=s.n; t++) {
X
X opr_qZ(t);
X
X qZ_mts += qZ;
X }
X}
X
X
Xrealmat QZ_mts;
X
Xvoid opr_QZ_mts()
X{
X QZ_mts.resize(s.M*s.K, s.p, REAL_ZERO );
X
X for (INTEGER t=1; t<=s.n; t++) {
X
X opr_QZ(t);
X
X QZ_mts += QZ;
X }
X}
X
X
Xrealmat ZZ;
X
Xvoid opr_ZZ()
X{
X ZZ.resize(s.K, s.K, REAL_ZERO );
X
X for (INTEGER t=1; t<=s.n; t++) {
X
X m.Z(t);
X
X ZZ += m.Z_t * T(m.Z_t);
X }
X}
X
X
Xvoid TSLS_obj()
X{
X INTEGER alpha, i, ii;
X INTEGER beta, j, jj;
X
X realmat varinv = invpsd(s.var,s.eps);
X
X opr_qZ_mts();
X opr_ZZ();
X
X realmat ZZinv = invpsd(ZZ,s.eps);
X s.obj = REAL_ZERO ;
X
X for (alpha=1; alpha<=s.M; alpha++) {
X for (i=1; i<=s.K; i++) {
X for (beta=1; beta<=s.M; beta++) {
X for (j=1; j<=s.K; j++) {
X
X ii = s.K*(alpha-1) + i;
X jj = s.K*(beta-1) + j;
X
X s.obj += qZ_mts[ii] * qZ_mts[jj]
X * varinv.elem(alpha,beta) * ZZinv.elem(i,j);
X }
X }
X }
X }
X}
X
X
Xvoid TSLS_mgn()
X{
X INTEGER alpha, i, ii;
X INTEGER beta, j, jj;
X INTEGER k,l;
X
X realmat varinv = invpsd(s.var,s.eps);
X
X opr_qZ_mts();
X opr_QZ_mts();
X opr_ZZ();
X
X realmat ZZinv = invpsd(ZZ,s.eps);
X
X s.obj = REAL_ZERO ;
X s.V.resize(s.p, s.p, REAL_ZERO );
X s.D.resize(s.p, 1, REAL_ZERO );
X
X for (alpha=1; alpha<=s.M; alpha++) {
X for (i=1; i<=s.K; i++) {
X for (beta=1; beta<=s.M; beta++) {
X for (j=1; j<=s.K; j++) {
X
X ii = s.K*(alpha-1) + i;
X jj = s.K*(beta-1) + j;
X
X s.obj += qZ_mts[ii] * qZ_mts[jj]
X * varinv.elem(alpha,beta) * ZZinv.elem(i,j);
X
X for (k=1; k<=s.p; k++) {
X
X s.D[k] += QZ_mts.elem(ii,k) * qZ_mts[jj]
X * varinv.elem(alpha,beta) * ZZinv.elem(i,j);
X
X for (l=1; l<=s.p; l++) {
X
X s.V.elem(k,l) += QZ_mts.elem(ii,k) * QZ_mts.elem(jj,l)
X * varinv.elem(alpha,beta) * ZZinv.elem(i,j);
X
X }
X }
X }
X }
X }
X }
X
X s.V = invpsd(s.V,s.eps);
X
X s.D = - s.V * s.D;
X}
X
X
Xvoid TSLS_V()
X{
X s.rank = 0;
X for (INTEGER i=1; i<=s.p; i++)
X if (s.V.elem(i,i) > REAL_ZERO ) s.rank++;
X};
X
X
Xvoid GMM_obj()
X{
X realmat sse;
X realmat varinv = invpsd(s.var,s.eps);
X
X opr_qZ_mts();
X
X sse = T(qZ_mts) * varinv * qZ_mts ;
X
X s.obj = sse[1];
X}
X
X
Xvoid GMM_mgn()
X{
X realmat sse;
X realmat varinv = invpsd(s.var,s.eps);
X
X opr_qZ_mts();
X opr_QZ_mts();
X
X sse = T(qZ_mts) * varinv * qZ_mts ;
X
X s.obj = sse[1];
X
X s.D = T(QZ_mts) * varinv * qZ_mts ;
X
X s.V = T(QZ_mts) * varinv * QZ_mts ;
X
X s.V = invpsd(s.V,s.eps);
X
X s.D = - s.V * s.D;
X}
X
X
Xvoid GMM_var(int var_loop)
X{
X INTEGER l = s.M * s.K;
X INTEGER t;
X INTEGER tau;
X INTEGER alpha, i, j, ii, jj;
X REAL x,weight;
X
X strcpy(s.df,"uncorrected");
X
X if (rows(s.var)==0 || cols(s.var)==0) {
X
X s.var.resize(l,l,REAL_ZERO );
X
X opr_ZZ();
X
X for (alpha=1; alpha<=s.M; alpha++) {
X for (i=1; i<=s.K; i++) {
X for (j=1; j<=s.K; j++) {
X
X ii = s.K*(alpha-1) + i;
X jj = s.K*(alpha-1) + j;
X
X s.var.elem(ii,jj) = ZZ.elem(i,j);
X
X }
X }
X }
X }
X if (var_loop == 0) return;
X
X realmat I(l,l,REAL_ZERO);
X realmat I_tau(l,l);
X realmat qZ_lag(l,1);
X
X for (tau=0; tau<=s.MA; tau++) {
X
X for (i=1; i<=l; i++)
X for (j=1; j<=l; j++)
X I_tau.elem(i,j)=REAL_ZERO ;
X
X for (t=tau+1; t<=s.n; t++) {
X
X opr_qZ(t-tau);
X qZ_lag = qZ;
X opr_qZ(t);
X
X I_tau += qZ * T(qZ_lag);
X }
X
X if (strcmp(s.weights,"Parzen") == 0 && tau > 0) {
X x = tau/(REAL)s.MA;
X if ( x < 0.5 )
X weight = 1.0 - 6.0*pow(x,2) + 6.0*pow(x,3);
X else
X weight = 2.0*pow((1.0 - x),3);
X
X }
X else {
X weight = 1.0;
X }
X
X I = I + weight*I_tau;
X if (tau > 0) {
X I = I + weight*T(I_tau);
X }
X }
X
X s.var = I;
X
X return;
X}
X
X
Xvoid GMM_V()
X{
X s.rank = 0;
X for (INTEGER i=1; i<=s.p; i++)
X if (s.V.elem(i,i) > REAL_ZERO ) s.rank++;
X}
END_OF_FILE
if test 9187 -ne `wc -c <'nlmdl.opr'`; then
echo shar: \"'nlmdl.opr'\" unpacked with wrong size!
fi
# end of 'nlmdl.opr'
fi
echo shar: End of archive 3 \(of 6\).
cp /dev/null ark3isdone
MISSING=""
for I in 1 2 3 4 5 6 ; do
if test ! -f ark${I}isdone ; then
MISSING="${MISSING} ${I}"
fi
done
if test "${MISSING}" = "" ; then
echo You have unpacked all 6 archives.
rm -f ark[1-9]isdone
else
echo You still need to unpack the following archives:
echo " " ${MISSING}
fi
## End of shell archive.
exit 0
exit 0 # Just in case...
--
Kent Landfield INTERNET: kent at sparky.IMD.Sterling.COM
Sterling Software, IMD UUCP: uunet!sparky!kent
Phone: (402) 291-8300 FAX: (402) 291-4362
Please send comp.sources.misc-related mail to kent at uunet.uu.net.
More information about the Comp.sources.misc
mailing list