# Inverse and backslash \\ division

**URL:** https://discourse.julialang.org/t/inverse-and-backslash-division/44680
**Category:** Statistics
**Tags:** question, optimization
**Created:** [August 10, 2020, 1:19pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680 "2020-08-10T13:19:31Z")
**Posts on this page:** 6
**Page:** 1

<div class="post-metadata">

### Author: ![luboshanus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/luboshanus/32/12163_2.png) [@luboshanus](https://discourse.julialang.org/u/luboshanus)
#### Post date: [August 10, 2020, 1:19pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/1 "2020-08-10T13:19:31Z")

</div>

Hi,

could anyone please explain to me how to use `inv()` and `\` properly to have valid results. Alse when checking with `GLM` package basic OLS gives different results and sometimes it is way too far. It might be caused by the data but strangely.

Thanks a lot.

```julia
x = randn(100, 3)
y = randn(100)

```

```julia
julia> (transpose(x) * x) \ (transpose(x) *y) .== inv(transpose(x) * x) * (transpose(x) *y)
3-element BitArray{1}:
1
0
1

julia> (transpose(x) * x) \ (transpose(x) *y)
3-element Array{Float64,1}:
-0.13268579177064818
-0.0009498174455286444
0.06340639348585342

julia> inv(transpose(x) * x) * (transpose(x) *y)
3-element Array{Float64,1}:
-0.13268579177064818
-0.0009498174455286456
0.06340639348585342

julia> using GLM

julia> (transpose(x) * x) \ (transpose(x) *y)
3-element Array{Float64,1}:
-0.13268579177064818
-0.0009498174455286444
0.06340639348585342

julia> GLM.coef(lm(x, y))
3-element Array{Float64,1}:
-0.1326857917706482
-0.0009498174455286441
0.06340639348585342

julia> (transpose(x) * x) \ (transpose(x) *y) .== GLM.coef(lm(x, y))
3-element BitArray{1}:
0
0
1

```

---

<div class="post-metadata">

### Author: ![dilumaluthge](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dilumaluthge/32/29283_2.png) [@dilumaluthge](https://discourse.julialang.org/u/dilumaluthge)
#### Post date: [August 10, 2020, 1:42pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/2 "2020-08-10T13:42:47Z")

</div>

Instead of checking for exact equality with `==`, I would recommend using approximate equality. See `isapprox`.

You should be able to just replace the `==` with `≈`.

See: [Mathematics · The Julia Language](https://docs.julialang.org/en/v1/base/math/#Base.isapprox)

---

<div class="post-metadata">

### Author: ![luboshanus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/luboshanus/32/12163_2.png) [@luboshanus](https://discourse.julialang.org/u/luboshanus)
#### Post date: [August 10, 2020, 1:51pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/3 "2020-08-10T13:51:43Z")

</div>

Thanks for an answer, but in my case it is more about the results, which, for data I have, lead to exploding results in OLS because of definiteness of `inv(X^T X)` matrix etc.

Then the difference is super significant not only on 6th+ decimal points but on 2nd, for example. And values really differ when using `inv` or `\`.

Any ideas? Thanks.

---

<div class="post-metadata">

### Author: ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [August 10, 2020, 2:28pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/4 "2020-08-10T14:28:28Z")

</div>

> [@luboshanus](#):
>
> Thanks for an answer, but in my case it is more about the results, which, for data I have, lead to exploding results in OLS because of definiteness of `inv(X^T X)` matrix etc.

It sounds like you are solving a least-squares problem by the “normal-equations” approach, i.e. you are computing z = (X^T X)^{-1} X^T y to minimize the least-square error \Vert X z - y \Vert\_2. If you are getting “exploding” results, it is presumably because X^T X is badly conditioned (nearly singular), in which case you have to be careful about how you solve the problem — your answer is very sensitive to roundoff errors and other errors (e.g. noise in the data).

If X^T X is badly conditioned, then `inv(X'X) * (X'y)` could indeed give very different answers than `(X'X) \ (X'y)` because of differences in roundoff errors — but the correct interpretation is that _both_ approaches are giving you an inaccurate result.

A better approach would be to use `X \ y`, which is equivalent to `qr(X, Val(true)) \ y` — this uses a pivoted QR factorization to solve the least-squares problem, which is much more accurate than the normal-equations approach for badly conditioned X.

If you need even more control over regularization for a badly conditioned problem, you could use a pseudo-inverse (`pinv`) approach. See also my explanation in another thread: [Efficient way of doing linear regression - #33 by stevengj](https://discourse.julialang.org/t/efficient-way-of-doing-linear-regression/31232/33)

---

<div class="post-metadata">

### Author: ![luboshanus](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/luboshanus/32/12163_2.png) [@luboshanus](https://discourse.julialang.org/u/luboshanus)
#### Post date: [August 10, 2020, 2:58pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/5 "2020-08-10T14:58:08Z")

</div>

Thanks. This I do understand. And will try to use QR factorization.

Also, I will add some information/comment to the problem I’m looking at. It is an analysis where I many perform OLS on expanding window over data. And the thing I have realised is that on 1:1000 observations OLS is ok, and when window is 1:1001, the coefficients explode. This is strange since the window is expanding and 1 observation pushes the analysis into unstable solution and when I move one observation further, 1:1002, it is stable again. So the matrix `XTX` is probably on the edge of well/bad conditioning, or I do not have any better explanation. Does this make sense in accordance to what you’ve written?

---

<div class="post-metadata">

### Author: ![chakravala](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chakravala/32/6832_2.png) [@chakravala](https://discourse.julialang.org/u/chakravala)
#### Post date: [August 14, 2020, 8:18pm UTC](https://discourse.julialang.org/t/inverse-and-backslash-division/44680/6 "2020-08-14T20:18:06Z")

</div>

I have come up with a more efficient method than the SaticArrays method for an underdetermined system:

> [@Grassmann.jl A\\b 3x faster than Julia's StaticArrays.jl](https://discourse.julialang.org/t/grassmann-jl-a-b-3x-faster-than-julias-staticarrays-jl/41451/35):
>
> This was because of a special explicit case for dimensions 1 and 2 and 3, which are now accounted for in Grassmann.jl also. Also, support has been added for Moore-Penrose inverses for underdetermined and overdetermined linear systems. For underdetermined cases, the exterior product algorithm works ~20x faster than the SMatrix algorithm, and it is numerically stable. For overdetermined equations, the method used is based on the traditional normal equations, and this is prone to more numerical i…

For underdetermined systems, my new method is about ~20x faster than the `SMatrix` and it is numerically stable. However, for overdetermined cases, I am falling back on the traditional normal equations, which are less numerically stable.
