You should never, ever do OLS by calculating (X'X)^{-1} X'Y: it can be
ill-conditioned even for smallish matrices in toy problems.Use QR factorization (or SVD in some special cases, especially if you need to regularize). Any numerical text should have the details, but fortunately you don't even need to program the QR solution in Julia, since X \ Y does the Right Thing™. See http://julia.readthedocs.org/en/latest/stdlib/linalg/ Best, Tamas On Fri, Jan 02 2015, Gabriel Mihalache <[email protected]> wrote: > AR is just OLS, so you should be able to do > > (X' * X) \ X' * Y > > with the right X and Y matrices. > > Note the difference between * and .* OLS is a matrix expression, not a > element-wise operation. > > Start by making sure that the X and Y matrices look like you expect them > to, shape-wise and content-wise. > > Regarding the TimeSeries package, I can't test my hunch right now, but I > suspect it's because 0 and 0.0 are different in Julia. One's an Integer and > the other's a Float. Same with 1 and 1.0 My guess is that the TimeSeries > function needs Float arguments and you pass it Integers. > > On Thursday, January 1, 2015 7:58:05 PM UTC-5, jspark wrote: >> >> Hi, >> >> Now I become a two weeks old for Julia and still struggling with Julia. >> >> I am porting *AR1 function* below from a *cran R* but for some reason it >> gives me an error (SingularException(3)) when I convert *solve() to get >> Alpha* into *inv()*. However, it works with R and Matlab. >> >> y=float(r,1]) # r is data vector from the attached DJ.csv file >> >> T = length(y)-1; >> >> Y = y[2:T]'; >> >> Y1= y[1:(T-1)]'; >> >> Alpha= inv(Y1'.*Y).*(Y1').*Y; >> >> >> *ERROR : SingularException(3)* >> >> As an alternative, I am trying to get AR1 coefficient from *arima() *of >> the *TimeModels* package but it's not working either. >> >> using TimeModels >> >> zz=rand(500,1); #500*1 Array(Float64,2) >> arima(zz,1,0,0); >> >> *ERROR : 'arima' has no method matching arima (::Array{Float64,2), >> ::Int64, ::Int64, ::Int64)* >> >> Can someone tell me what is the behind story of the error? >> >> Many thanks and Happy New Year to all!! >> >> >> AR1 <-function(x){T <- length(x) -1Y <- x[2:T]Y_ <- x[1:(T-1)] ALPHA <- >> solve(t(Y_) %*% Y_ ) %*% t(Y_) %*% YRE <- Y - Y_ %*% ALPHASIGMAS <- sum(RE >> ^2) / (T-1)STDA <- sqrt( SIGMAS * solve(t(Y_) %*% Y_ ))return(list(ALPHA= >> ALPHA,STDA=STDA))} Source : ttps:// >> github.com/cran/vrtest/blob/master/R/AR1.R >> >>
