Mercurial > octave-nkf
annotate scripts/statistics/base/gls.m @ 20266:83792dd9bcc1
Use in-place operators in m-files where possible.
* scripts/audio/@audioplayer/set.m, scripts/audio/@audiorecorder/set.m,
scripts/audio/mu2lin.m, scripts/elfun/cosd.m, scripts/general/del2.m,
scripts/general/profexplore.m, scripts/general/quadl.m, scripts/general/rat.m,
scripts/general/rotdim.m, scripts/help/get_first_help_sentence.m,
scripts/help/private/__strip_html_tags__.m, scripts/image/cubehelix.m,
scripts/io/textread.m, scripts/linear-algebra/duplication_matrix.m,
scripts/linear-algebra/housh.m, scripts/linear-algebra/krylov.m,
scripts/linear-algebra/logm.m, scripts/linear-algebra/normest.m,
scripts/linear-algebra/onenormest.m, scripts/optimization/fminsearch.m,
scripts/optimization/lsqnonneg.m, scripts/optimization/qp.m,
scripts/plot/appearance/annotation.m, scripts/plot/appearance/axis.m,
scripts/plot/appearance/legend.m, scripts/plot/appearance/specular.m,
scripts/plot/draw/colorbar.m, scripts/plot/draw/hist.m,
scripts/plot/draw/plotmatrix.m, scripts/plot/draw/private/__stem__.m,
scripts/plot/util/__actual_axis_position__.m,
scripts/plot/util/__gnuplot_drawnow__.m, scripts/plot/util/findobj.m,
scripts/plot/util/print.m, scripts/plot/util/private/__go_draw_axes__.m,
scripts/plot/util/private/__print_parse_opts__.m, scripts/plot/util/rotate.m,
scripts/polynomial/pchip.m, scripts/polynomial/polyaffine.m,
scripts/polynomial/polyder.m, scripts/polynomial/private/__splinefit__.m,
scripts/polynomial/residue.m, scripts/signal/arch_fit.m,
scripts/signal/arch_rnd.m, scripts/signal/bartlett.m,
scripts/signal/blackman.m, scripts/signal/freqz.m, scripts/signal/hamming.m,
scripts/signal/hanning.m, scripts/signal/spectral_adf.m,
scripts/signal/spectral_xdf.m, scripts/signal/stft.m,
scripts/sparse/bicgstab.m, scripts/sparse/cgs.m,
scripts/sparse/private/__sprand_impl__.m, scripts/sparse/qmr.m,
scripts/sparse/sprandsym.m, scripts/sparse/svds.m, scripts/specfun/legendre.m,
scripts/special-matrix/gallery.m, scripts/statistics/base/gls.m,
scripts/statistics/models/logistic_regression.m,
scripts/statistics/tests/kruskal_wallis_test.m,
scripts/statistics/tests/manova.m, scripts/statistics/tests/wilcoxon_test.m,
scripts/time/datevec.m:
Use in-place operators in m-files where possible.
| author | Rik <rik@octave.org> |
|---|---|
| date | Tue, 26 May 2015 21:07:42 -0700 |
| parents | b0f7ee81d974 |
| children |
| rev | line source |
|---|---|
|
19731
4197fc428c7d
maint: Update copyright notices for 2015.
John W. Eaton <jwe@octave.org>
parents:
17744
diff
changeset
|
1 ## Copyright (C) 1996-2015 John W. Eaton |
| 3200 | 2 ## |
| 3 ## This file is part of Octave. | |
| 4 ## | |
| 5 ## Octave is free software; you can redistribute it and/or modify it | |
| 6 ## under the terms of the GNU General Public License as published by | |
| 7016 | 7 ## the Free Software Foundation; either version 3 of the License, or (at |
| 8 ## your option) any later version. | |
| 3200 | 9 ## |
| 10 ## Octave is distributed in the hope that it will be useful, but | |
| 11 ## WITHOUT ANY WARRANTY; without even the implied warranty of | |
| 12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | |
| 13 ## General Public License for more details. | |
| 14 ## | |
| 15 ## You should have received a copy of the GNU General Public License | |
| 7016 | 16 ## along with Octave; see the file COPYING. If not, see |
| 17 ## <http://www.gnu.org/licenses/>. | |
| 3200 | 18 |
| 3381 | 19 ## -*- texinfo -*- |
| 3368 | 20 ## @deftypefn {Function File} {[@var{beta}, @var{v}, @var{r}] =} gls (@var{y}, @var{x}, @var{o}) |
|
20209
d9341b422488
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19867
diff
changeset
|
21 ## Generalized least squares model. |
|
d9341b422488
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19867
diff
changeset
|
22 ## |
|
d9341b422488
doc: Update more docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents:
19867
diff
changeset
|
23 ## Perform a generalized least squares estimation for the multivariate model |
| 3368 | 24 ## @tex |
| 25 ## $y = x b + e$ | |
| 26 ## with $\bar{e} = 0$ and cov(vec($e$)) = $(s^2)o$, | |
| 27 ## @end tex | |
|
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7017
diff
changeset
|
28 ## @ifnottex |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
29 ## @w{@math{y = x*b + e}} with @math{mean (e) = 0} and |
| 3499 | 30 ## @math{cov (vec (e)) = (s^2) o}, |
|
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7017
diff
changeset
|
31 ## @end ifnottex |
|
20229
b0f7ee81d974
doc: Remove extra spaces at start of docstring which show up in Info format.
Rik <rik@octave.org>
parents:
20209
diff
changeset
|
32 ## where |
| 3426 | 33 ## @tex |
| 3368 | 34 ## $y$ is a $t \times p$ matrix, $x$ is a $t \times k$ matrix, $b$ is a $k |
| 35 ## \times p$ matrix, $e$ is a $t \times p$ matrix, and $o$ is a $tp \times | |
| 36 ## tp$ matrix. | |
| 37 ## @end tex | |
|
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7017
diff
changeset
|
38 ## @ifnottex |
| 3499 | 39 ## @math{y} is a @math{t} by @math{p} matrix, @math{x} is a @math{t} by |
| 40 ## @math{k} matrix, @math{b} is a @math{k} by @math{p} matrix, @math{e} | |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
41 ## is a @math{t} by @math{p} matrix, and @math{o} is a @math{t*p} by |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
42 ## @math{t*p} matrix. |
|
8517
81d6ab3ac93c
Allow documentation tobe built for other formats than tex and info
sh@sh-laptop
parents:
7017
diff
changeset
|
43 ## @end ifnottex |
| 3426 | 44 ## |
| 3368 | 45 ## @noindent |
| 3499 | 46 ## Each row of @var{y} and @var{x} is an observation and each column a |
| 47 ## variable. The return values @var{beta}, @var{v}, and @var{r} are | |
| 48 ## defined as follows. | |
| 3426 | 49 ## |
| 3368 | 50 ## @table @var |
| 51 ## @item beta | |
| 3499 | 52 ## The GLS estimator for @math{b}. |
| 3426 | 53 ## |
| 3368 | 54 ## @item v |
| 3499 | 55 ## The GLS estimator for @math{s^2}. |
| 3426 | 56 ## |
| 3368 | 57 ## @item r |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
58 ## The matrix of GLS residuals, @math{r = y - x*beta}. |
| 3368 | 59 ## @end table |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
60 ## @seealso{ols} |
| 3368 | 61 ## @end deftypefn |
| 3200 | 62 |
| 63 ## Author: Teresa Twaroch <twaroch@ci.tuwien.ac.at> | |
| 64 ## Created: May 1993 | |
| 65 ## Adapted-By: jwe | |
| 66 | |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
67 function [beta, v, r] = gls (y, x, o) |
| 3200 | 68 |
| 69 if (nargin != 3) | |
| 6046 | 70 print_usage (); |
| 3200 | 71 endif |
| 72 | |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
73 if (! (isnumeric (x) && isnumeric (y) && isnumeric (o))) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
74 error ("gls: X, Y, and O must be numeric matrices or vectors"); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
75 endif |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
76 |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
77 if (ndims (x) != 2 || ndims (y) != 2 || ndims (o) != 2) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
78 error ("gls: X, Y and O must be 2-D matrices or vectors"); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
79 endif |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
80 |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
81 [rx, cx] = size (x); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
82 [ry, cy] = size (y); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
83 [ro, co] = size (o); |
| 3200 | 84 if (rx != ry) |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
85 error ("gls: number of rows of X and Y must be equal"); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
86 endif |
|
19867
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19731
diff
changeset
|
87 if (! issquare (o) || ro != ry*cy) |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
88 error ("gls: matrix O must be square matrix with rows = rows (Y) * cols (Y)"); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
89 endif |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
90 |
|
12656
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
91 if (isinteger (x)) |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
92 x = double (x); |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
93 endif |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
94 if (isinteger (y)) |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
95 y = double (y); |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
96 endif |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
97 if (isinteger (o)) |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
98 o = double (o); |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
99 endif |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
100 |
|
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
101 ## Start of algorithm |
|
20266
83792dd9bcc1
Use in-place operators in m-files where possible.
Rik <rik@octave.org>
parents:
20229
diff
changeset
|
102 o ^= -1/2; |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
103 z = kron (eye (cy), x); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
104 z = o * z; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
105 y1 = o * reshape (y, ry*cy, 1); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
106 u = z' * z; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
107 r = rank (u); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
108 |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
109 if (r == cx*cy) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
110 b = inv (u) * z' * y1; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
111 else |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
112 b = pinv (z) * y1; |
| 3200 | 113 endif |
| 114 | |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
115 beta = reshape (b, cx, cy); |
| 3200 | 116 |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
117 if (isargout (2) || isargout (3)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
118 r = y - x * beta; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
119 if (isargout (2)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
120 v = (reshape (r, ry*cy, 1))' * (o^2) * reshape (r, ry*cy, 1) / (rx*cy - r); |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
121 endif |
| 3200 | 122 endif |
| 123 | |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
124 endfunction |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
125 |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
126 |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
127 %!test |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
128 %! x = [1:5]'; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
129 %! y = 3*x + 2; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
130 %! x = [x, ones(5,1)]; |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
131 %! o = diag (ones (5,1)); |
|
12656
6b2f14af2360
Overhaul functions in statistics/base directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
132 %! assert (gls (y,x,o), [3; 2], 50*eps); |
| 3200 | 133 |
|
19867
9fc020886ae9
maint: Clean up m-files to follow Octave coding conventions.
Rik <rik@octave.org>
parents:
19731
diff
changeset
|
134 ## Test input validation |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
135 %!error gls () |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
136 %!error gls (1) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
137 %!error gls (1, 2) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
138 %!error gls (1, 2, 3, 4) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
139 %!error gls ([true, true], [1, 2], ones (2)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
140 %!error gls ([1, 2], [true, true], ones (2)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
141 %!error gls ([1, 2], [1, 2], true (2)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
142 %!error gls (ones (2,2,2), ones (2,2), ones (4,4)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
143 %!error gls (ones (2,2), ones (2,2,2), ones (4,4)) |
|
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
144 %!error gls (ones (2,2), ones (2,2), ones (4,4,4)) |
|
14363
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
145 %!error gls (ones (1,2), ones (2,2), ones (2,2)) |
|
f3d52523cde1
Use Octave coding conventions in all m-file %!test blocks
Rik <octave@nomad.inbox5.com>
parents:
14138
diff
changeset
|
146 %!error gls (ones (2,2), ones (2,2), ones (2,2)) |
|
11436
e151e23f73bc
Overhaul base statistics functions and documentation of same.
Rik <octave@nomad.inbox5.com>
parents:
9211
diff
changeset
|
147 |