# Strange error when performing linear regression

**URL:** <https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303>\
**Category:** General Usage\
**Created:** [December 18, 2021, 7:35pm UTC](https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303 "2021-12-18T19:35:08Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![rkube](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rkube/32/211198_2.png) [@rkube](https://discourse.julialang.org/u/rkube)\
**Post date:** [December 18, 2021, 7:35pm UTC](https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303/1 "2021-12-18T19:35:08Z")

</div>

What is going on here?

```julia
julia> x = log.(Nz_list[1:end])
5-element Vector{Float64}:
 4.1588830833596715
 4.852030263919617
 5.545177444479562
 6.238324625039508
 6.931471805599453

julia> 

julia> y = log.(dist_list[1:end])
5-element Vector{Float64}:
 -4.074961466558573
 -5.1239528564742365
 -6.165969558863518
 -7.206262977477147
 -8.246126826893256

julia> (x' * x) \ (x' * y)
-1.1233909121590568

julia> x = log.(Nz_list)
1×5 Matrix{Float64}:
 4.15888 4.85203 5.54518 6.23832 6.93147

julia> y = log.(dist_list)
1×5 Matrix{Float64}:
 -4.07496 -5.12395 -6.16597 -7.20626 -8.24613

julia> (x' * x) \ (x' * y)
ERROR: SingularException(4)
Stacktrace:
 [1] checknonsingular
   @ /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:19 [inlined]
 [2] checknonsingular
   @ /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/factorization.jl:21 [inlined]
 [3] #lu!#136
   @ /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:85 [inlined]
 [4] #lu#140
   @ /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:273 [inlined]
 [5] lu (repeats 2 times)
   @ /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/lu.jl:272 [inlined]
 [6] \(A::Matrix{Float64}, B::Matrix{Float64})
   @ LinearAlgebra /builddir/build/BUILD/julia-1.6.4/build/usr/share/julia/stdlib/v1.6/LinearAlgebra/src/generic.jl:1136
 [7] top-level scope
   @ REPL[89]:1

```

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [December 18, 2021, 7:56pm UTC](https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303/2 "2021-12-18T19:56:03Z")

</div>

Don’t use the normal equations. Just use `x \ y` which will use a QR factorization and be more accurate.

---

<div class="post-metadata">

**Author:** ![oheil](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oheil/32/220745_2.png) [@oheil](https://discourse.julialang.org/u/oheil)\
**Post date:** [December 18, 2021, 7:56pm UTC](https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303/3 "2021-12-18T19:56:11Z")

</div>

Not an expert of Linear Algebra here, but if I remember correct:  
`(x' * x)` is a scalar in the case of x being a vector, but a matrix in the case of x being a matrix:

```julia
julia> x=reshape([4.15888, 4.85203, 5.54518, 6.23832, 6.93147], 1, 5 )
1×5 Matrix{Float64}:
 4.15888 4.85203 5.54518 6.23832 6.93147

julia> x' * x
5×5 Matrix{Float64}:
 17.2963 20.179 23.0617 25.9444 28.8272
 20.179 23.5422 26.9054 30.2685 33.6317
 23.0617 26.9054 30.749 34.5926 38.4362
 25.9444 30.2685 34.5926 38.9166 43.2407
 28.8272 33.6317 38.4362 43.2407 48.0453

```

versus:

```julia
julia> x=[4.15888, 4.85203, 5.54518, 6.23832, 6.93147]
5-element Vector{Float64}:
 4.15888
 4.85203
 5.54518
 6.23832
 6.93147

julia> x' * x
158.549411991

```

Is this good enough for an explanation? Well, a real expert is replying… lets wait 😉

---

<div class="post-metadata">

**Author:** ![John\_Gibson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/john_gibson/32/5321_2.png) [@John\_Gibson](https://discourse.julialang.org/u/John_Gibson)\
**Post date:** [December 18, 2021, 10:01pm UTC](https://discourse.julialang.org/t/strange-error-when-performing-linear-regression/73303/4 "2021-12-18T22:01:39Z")

</div>

Exactly. In the first case where `x` and `y` are vectors, `x' * x` and `y' *y` are scalar valued inner products, and `(x' * x) \ (x' * y)` is scalar division.

```julia
julia> x = [4.1588830833596715, 4.852030263919617, 5.545177444479562, 6.238324625039508, 6.931471805599453]
5-element Vector{Float64}:
 4.1588830833596715
 4.852030263919617
 5.545177444479562
 6.238324625039508
 6.931471805599453

julia> y = [-4.074961466558573, -5.1239528564742365, -6.165969558863518, -7.206262977477147, -8.246126826893256]
5-element Vector{Float64}:
 -4.074961466558573
 -5.1239528564742365
 -6.165969558863518
 -7.206262977477147
 -8.246126826893256

julia> x' * x
158.54949459300647

julia> x' * y
-178.113061353195

julia> (x' * x) \ (x' * y)
-1.1233909121590568

```

In the second case where `x` is a 1 x 5 matrix, `x' * x` is a 5x1 x 1x5 outer product matrix, which is necessarily rank-deficient. Let `A = x' * x`. Then the jth column of A is the vector `x'` times the scalar `x[j]`. So all the columns of `A` are linearly dependent, and `A` has rank 1.

```julia
julia> x = x'
1×5 adjoint(::Vector{Float64}) with eltype Float64:
 4.15888 4.85203 5.54518 6.23832 6.93147

julia> y = y'
1×5 adjoint(::Vector{Float64}) with eltype Float64:
 -4.07496 -5.12395 -6.16597 -7.20626 -8.24613

julia> A = x' * x
5×5 Matrix{Float64}:
 17.2963 20.179 23.0617 25.9445 28.8272
 20.179 23.5422 26.9054 30.2685 33.6317
 23.0617 26.9054 30.749 34.5926 38.4362
 25.9445 30.2685 34.5926 38.9167 43.2408
 28.8272 33.6317 38.4362 43.2408 48.0453

julia> A[:,1]
5-element Vector{Float64}:
 17.296308501055247
 20.179026584564458
 23.061744668073665
 25.944462751582876
 28.827180835092083

julia> x' * x[1]
5-element Vector{Float64}:
 17.296308501055247
 20.179026584564458
 23.061744668073665
 25.944462751582876
 28.827180835092083

julia> using LinearAlgebra

julia> rank(A)
1

```