COMPGROUPS.NET | Search | Post Question | Groups | Stream | About | Register

• Email
• Follow

```Hi everybody,

-I'm trying to program the ordinary least squares method, and I must
solve this matrix system:

e = h - A* x
x is the unknown.

After the resolution of the system and when I try the verification, x
doesn't resolve it!
Why?

The program:
% ------------------------------------------
% ------------------------------------------
clc;

M = 1; % arbitrary
N = 5; % arbitrary
P = 5; % arbitrary but P >= Q
Q = 4; % arbitrary but P >= Q

e = zeros (M* (N+ P+ 1), 1); % arbitrary = 0
h = ones (M* (N+ P+ 1), 1); % arbitrary = ones
A = eye (M* (N+ P+ 1), P+ M* (Q+ 1)); % arbitrary = eye

x = ((A' * A)^(-1)) * (A' * (h- e)),
% the size of x is (P+ M* (Q+ 1), 1);
% (A' * A) to transform into a square matrix for inversion

e == (h - (A * x)),
% Normally after this, the result is 1
% Why it doesn't!?

% ------------------------------------------
% ------------------------------------------

```
 0
Reply zaki_nabil (105) 3/16/2005 8:13:09 PM

See related articles to this posting

```In article <eeff36b.-1@webx.raydaftYaTP>, "Zaki N."
<zaki_nabil@caramail.com> wrote:

> Hi everybody,
>
> -I'm trying to program the ordinary least squares method, and I must
> solve this matrix system:
>
> e = h - A* x
> x is the unknown.
>
> After the resolution of the system and when I try the verification, x
> doesn't resolve it!
> Why?
>
> The program:
> % ------------------------------------------
> clc;
>
> M = 1; % arbitrary
> N = 5; % arbitrary
> P = 5; % arbitrary but P >= Q
> Q = 4; % arbitrary but P >= Q
>
> e = zeros (M* (N+ P+ 1), 1); % arbitrary = 0
> h = ones (M* (N+ P+ 1), 1); % arbitrary = ones
> A = eye (M* (N+ P+ 1), P+ M* (Q+ 1)); % arbitrary = eye
>
> x = ((A' * A)^(-1)) * (A' * (h- e)),
> % the size of x is (P+ M* (Q+ 1), 1);
> % (A' * A) to transform into a square matrix for inversion
>
> e == (h - (A * x)),
> % Normally after this, the result is 1
> % Why it doesn't!?
> % ------------------------------------------
>
-----------
Hello Zaki,

Two things are wrong here.  To begin with, you're requiring

e == (h - (A * x))

to be exactly true, which makes no allowance for round off errors.  You
need to change this to allow an error within a certain judiciously
selected tolerance.

However, the main problem is that, in effect, there are eleven linear
equalities here but only ten variables in x to adjust.  The idea in least
squares is to find the value of x which produces the smallest sum of
squares in the expression A*x-(h-e).  In the above circumstances it will,
in general, be impossible to make it identically zero.

(Remove "xyzzy" and ".invalid" to send me email.)
Roger Stafford
```
 0
Reply ellieandrogerxyzzy (4805) 3/16/2005 9:09:56 PM

```Thanks for your reply Roger.

-The problem isn't coming from the ordinary least squares method or
not.
We are in front of a problem of resolving a matrix system:

e = h - A* x

Where "x" is the unknown, and "h", "A", and "e" can be 'anything'
known.

After the resolution of the system and when I try the verification, x
doesn't resolve it!

Where is the mistake?
Please give me another method to resolve it.

-The program:

% ------------------------------------------
% ------------------------------------------
clc;

M = 1; % arbitrary
N = 5; % arbitrary
P = 5; % arbitrary but P >= Q
Q = 4; % arbitrary but P >= Q

e = zeros (M* (N+ P+ 1), 1); % arbitrary = 0
h = ones (M* (N+ P+ 1), 1); % arbitrary = ones
A = eye (M* (N+ P+ 1), P+ M* (Q+ 1)); % arbitrary = eye

x = ((A' * A)^(-1)) * (A' * (h- e)),
% the size of x is (P+ M* (Q+ 1), 1);
% (A' * A) to transform into a square matrix for inversion

e == (h - (A * x)),
% Normally after this, the result is 1
% Why it doesn't!?

% ------------------------------------------
% ------------------------------------------
```
 0
Reply zaki_nabil (105) 3/17/2005 7:40:37 AM

```"Zaki N." <zaki_nabil@caramail.com> wrote in message
news:eeff36b.1@webx.raydaftYaTP...
>
> -The problem isn't coming from the ordinary least squares method or
> not.
> We are in front of a problem of resolving a matrix system:
>
> e = h - A* x
>
> Where "x" is the unknown, and "h", "A", and "e" can be 'anything'
> known.
>
> After the resolution of the system and when I try the verification, x
> doesn't resolve it!

For there to be a solution x to the system A*x = h (i.e. the e vector is
exactly the zero vector) you would need to be able to write h as a linear
combination of the columns of A -- the coefficients in the linear
combination are the elements of x.  Note that no column of A contains a
nonzero value in the last element.  Therefore, no matter what coefficients
you choose, the linear combination of the columns of A will always have a
zero as the last element of the resulting vector.  Since h has a 1 as its
last element, we can't express h as a linear combination of the columns of A
and so the system A*x = h has no exact solution.

Since there is no solution x that satisfies A*x = h exactly, we try to solve
it in a least-squares sense -- minimize the norm of (A*x - h) [which is -e
in your terminology.]  In this case, the x that minimizes this norm is the
ones vector that the code below produces.  The residual (A*x - h) is:

[0; 0; 0; 0; 0; 0; 0; 0; 0; 0; -1]

The last element of the residual is always going to be -1, from the
explanation above.  Therefore to minimize the norm, we minimize the
magnitudes of the other elements of the residual -- and 0 is as small as the
magnitudes can get.

> Where is the mistake?

There is no mistake.

> Please give me another method to resolve it.

Just a suggestion:  don't use ^(-1) to invert your matrix.  If you need to
invert the matrix, use the INV function.  However, you don't need to invert
the matrix -- simply use the backslash operator (\) instead.  See HELP

--
Steve Lord
slord@mathworks.com

>
> -The program:
>
> % ------------------------------------------
> % ------------------------------------------
> clc;
>
> M = 1; % arbitrary
> N = 5; % arbitrary
> P = 5; % arbitrary but P >= Q
> Q = 4; % arbitrary but P >= Q
>
> e = zeros (M* (N+ P+ 1), 1); % arbitrary = 0
> h = ones (M* (N+ P+ 1), 1); % arbitrary = ones
> A = eye (M* (N+ P+ 1), P+ M* (Q+ 1)); % arbitrary = eye
>
> x = ((A' * A)^(-1)) * (A' * (h- e)),
> % the size of x is (P+ M* (Q+ 1), 1);
> % (A' * A) to transform into a square matrix for inversion
>
> e == (h - (A * x)),
> % Normally after this, the result is 1
> % Why it doesn't!?
>
> % ------------------------------------------
> % ------------------------------------------

```
 0
Reply slord (13686) 3/17/2005 2:30:59 PM

```Thank you!

Even if I didn't understood all your explanations!
```
 0
Reply zaki_nabil (105) 3/18/2005 9:31:36 AM

4 Replies
36 Views

Similar Articles

12/7/2013 5:13:18 PM
[PageSpeed]