Re: Large control loops
- To: mathgroup at smc.vnet.net
- Subject: [mg122355] Re: Large control loops
- From: Ray Koopman <koopman at sfu.ca>
- Date: Wed, 26 Oct 2011 17:38:58 -0400 (EDT)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- References: <j80qff$ael$1@smc.vnet.net> <j83aon$k7c$1@smc.vnet.net> <j862ol$5tl$1@smc.vnet.net>
On Oct 25, 3:24 am, Arthur <shuk... at gmail.com> wrote: > The vector in the first column is the common vector. It is the > independent variable in the subsequent regressions - I am > regressing first to find the correlation between it and each of > the other vectors. It's all numbers. > > OLS regressions. > > I am estimating the parameters through a secondary regression on > the residuals of the first and then some basic algebra. It's an > Ornstein-Uhlenbeck process, I am using the first procedure outlined > here > http://www.sitmo.com/article/calibrating-the-ornstein-uhlenbeck-model/ > > I am not so much worried about the time (although that too is a > concern) but the actual structure. Should it just be a long while > loop? par[r] expects r to be a vector of residuals. It uses van den Berg's first procedure and returns {a,b,sd,lambda,mu,sigma}. Since the residuals sum to zero, the formulas simplify a little. par[r_] := Block[{n = Length@r - 1, delta = 1 (* ?? *), Sx,Sy, Sxx,Syy, Sxy, a,b,sd}, {Sx , Sy } = -r[[{-1,1}]]; {Sxx, Syy} = n(#.#&@Take[r,{2,-2}] + {Sy,Sx}^2) - {Sx,Sy}^2 Sxy = n Most at r.Rest@r - Sx*Sy; { a = Sxy/Sxx, b = (Sy - a*Sx)/n, sd = Sqrt[(Syy - a*Sxy)/(n(n-2))], (* lambda = *) -Log[a]/delta, (* mu = *) b/(1-a), (* sigma = *) sd*Sqrt[-2 Log[a]/(delta(1-a^2))]} ] I would get all the residuals simultaneously, and then map par over them. {x,y} = {First@#,Last@#}&[Transpose@A - Mean@A]; p = par /@ (y + Transpose[{y.x/-x.x}].{x}]) Those two lines would be followed by a loop that scans p looking for the best parameters, and all that would be nested within a loop that puts different values into A. (Note: rows 1,2,... of p correspond to column 2,3,... of A.)