Mercurial > forge

changeset 10259:fde06f2c60b3 octave-forge

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
control-devel: arx: use Ljung's QR method
author paramaniac
date Tue, 15 May 2012 15:22:00 +0000
parents 14e70818755a
children 00bf15ad8b27
files extra/control-devel/inst/arx.m
diffstat 1 files changed, 10 insertions(+), 7 deletions(-) [+]
line wrap: on
line diff
--- a/extra/control-devel/inst/arx.m	Tue May 15 14:40:12 2012 +0000
+++ b/extra/control-devel/inst/arx.m	Tue May 15 15:22:00 2012 +0000
@@ -95,14 +95,17 @@
 
 
 function theta = __theta__ (phi, y, i, n)
-                                                
-  ## Theta = Phi \ Y(n+1:end, :);                   # naive formula
     
-  if (numel (phi) == 1)
-    ## single-experiment dataset
-    theta = __ls_svd__ (phi{1}, y{1}(n(i)+1:end, i));
-  else
-    ## multi-experiment dataset
+  if (numel (phi) == 1)                             # single-experiment dataset
+    A = horzcat (phi{1}, y{1}(n(i)+1:end, i));      # [Phi, Y]
+    R0 = triu (qr (A, 0));                          # 0 for economy-size R (without zero rows)
+    R1 = R0(1:end-1, 1:end-1);                      # R1 is triangular - can we exploit this in R1\R2
+    R2 = R0(1:end-1, end);
+    theta = __ls_svd__ (R1, R2);                    # R1 \ R2
+    
+    ## Theta = Phi \ Y(n+1:end, :);                 # naive formula
+    ## theta = __ls_svd__ (phi{1}, y{1}(n(i)+1:end, i));
+  else                                              # multi-experiment dataset
     ## TODO: find more sophisticated formula than
     ## Theta = (Phi1' Phi + Phi2' Phi2 + ...) \ (Phi1' Y1 + Phi2' Y2 + ...)