Matrix is singular, RCOND=NAN

I am using the following code in a simulation and get a warning message. Is there any way to correct it?

       lb = 0.00001*ones(2,1);
       ub = 0.99998*ones(2,1);

       options =optimset('Display','off','LargeScale','off');

       coef = lsqlin(bC,by,[],[],[],[], lb ,ub ,[],options);

Warning: Matrix is singular, close to singular or badly scaled. Results may be inaccurate. RCOND = NaN. > In qpsub at 359 In lsqlin at 336 In seasonalarch1 at 330

Thanks

Answer 1

[Matt J]

You are constructing C and y inadequately. You do not have enough linear equations, or they are not linearly independent enough, to uniquely identify the actual coef. Find additional and more informative equations and add them as new rows to C and y until cond(C) improves.

[dav]

All I have is the data set and I am dividing it to 5 vectors each containing 200 observation. if mod = 1 that data goes to yt1, if mod = 2 that data goes to yt2 and so on.. I dont understand what you mean by "Find more linear equations"

Also C and y contain data so Can I really add linear equations to those?? Could you explain it a little bit more please? may be an example

Thanks

[Matt J]

Here's a simpler example. The following equations are not linearly independent enough to give a unique solution for a and b

 a+b=0
 2*a+2*b=0

Like in your case, it is true that a solution very different from zero exists, e.g. a=-1,b=1. However, a=b=0 is also a solution.

By adding another equation, e.g. a+2*b=1, which is linearly independent of the initial ones, I can solve uniquely:

     >> C=[1 1; 2 2; 1 2]; y=[0;0;1];   coef=C\y; a=coef(1), b=coef(2)

    a =

       -1.0000

    b =

         1

All I have is the data set and I am dividing it to 5 vectors each containing 200 observation

It isn't enough. You must find more information about the coeffs from somewhere.

Also C and y contain data so Can I really add linear equations to those??

Yes, as my example showed.

[dav]

Thanks. However the only other information I am allowed to use is that the coefficients are greater than zero and the ar1 coef should be less than one.

I am using the first info here, even if I include ar1 coef <1, I get the same error!

[Walter Roberson]

Not every set of equations has a usable solution.

[Matt J]

I am using the first info here, even if I include ar1 coef <1, I get the same error!

That should come as no surprise. Again, if you look at my simplified example

 a+b=0
 2*a+2*b=0

this system still has infinite solutions even if you constrain 0<=a,b<=1.

If you really believe the data you have should be sufficiently informative, it is possible that you have some bug in your construction of C. For example, if x were accidentally made all zero, then all the rows of C would be identically [1,0] and of course no useful solution could come from that.

[Matt J]

In fact, the problem has to be occurring because x in the line

 C=[C,x]

is approximately constant (i.e., all x(i) are the same). Otherwise, C should be very well conditioned.

[dav]

I checked my C matrix and it was ok. I also checked the loops that create the x vectors and found nothing wrong.

When I change 0.00001 and 0.99998 to 0.01 and 0.98 I get the results I want. But, my worry is that I cannot change limits like that because the constraint is that the coefficients are between 0 and 1, not between 0.01 and 0.98.

[Matt J]

How did you check C? What was cond(C)?

Better yet, show us the 2x2 matrix C.'*C and the 2x1 vector C.'*y, so we can study the problem ourselves.

[dav]

I printed C and checked it to see if it has the problems you have mentioned earlier. I think you can change the above code a little bit to print 2x2 matrix C.'*C and the 2x1 vector C.'*y, if you want I can ajjust the above code to generate it and paste the new code here.

Thanks

[Matt J]

I don't have LSQLIN, but I can run rest of the code myself, and I do see now that C is okay. However, this means that I also cannot reproduce the problem. When I just apply an unconstrained solver

 coef = C\y

I find that coef is very different from zero and almost always satisfies the constraints, even without enforcing them specifically. What do you see when you do C\y?

In any case, you need to provide an example set of data which triggers the malfunction, so I can study what's going on without LSQLIN. For this, I don't need the entire C and y data matrices. Only C.'*C and C.'*y

[dav]

I figured out that this happens when the parameter estimates get closer to the boundary values I am using. sometimes LS estimates are still 0 and 1 even if I use a lower bound and an upper bound.

Theoretically, the estimates (if correct) should be greater than zero and less than oone. Can you please explain why lsqlin gives me 0 and 1 as estimate (sometimes- when i run it multiple times) even when I have a lb and an ub.

Thanks

[Matt J]

Still waiting for C.'*C and C.'*y

[dav]

I am sorry. here they are

>> C'*C

      ans =

        200.0000   34.0789
         34.0789   19.7842

    >> C'*y

    ans =

       30.3262
       12.7167

Thanks

[Matt J]

Can you please also give me the solution that LSQLIN is giving you for this data?

Also, what happens when you run as follows, instead

 coef = quadprog(C.'*C,-C.'*y,[],[],[],[],lb,ub);

[dav]

Since I had closed the program, I ran it again. The following results are for y = yt1 and C = [1 xt1]

alpha0 =

      0.5392

ar1 =

     1.0000e-05

>> C'*C

ans =

    200.0000   27.5925
     27.5925   18.4254

>> C'*y

ans =

    107.8483
     11.7406

As you can see ar coefficient is really small. this time it took the boundary value, but the actual value is 0.2 and sometimes it become equal to zero. Thanks

[dav]

Matt, could you also please let me know if I am using the fmincon code corresctly?

Tanks

   y = yt1;
        x= xt1;
        len = length(y);
        C = zeros(len,p);

        C(:,1) = 1; % The first column is for a0
        %for i = 1:p %Then we create shifted columns ( p in number ) for a 
         %  C(1+i:len,1+i) = x(1:len-i,1);     
        %end

        C = [C x];

         lb = [0.00001 0.00001]';
         ub = [0.99998 0.99998]';

        * fun=@(x) norm(C*x(:)-y)^2;*

         options =optimset('Display','off','LargeScale','off');

         *coef = fmincon(fun,y,[],[],[],[], lb ,ub,[],options);*

[Matt J]

As you can see ar coefficient is really small. this time it took the boundary value, but the actual value is 0.2 and sometimes it become equal to zero.

dav,

LSQLIN is giving you the correct values (for the problem data that you're providing it). Why do you think the correct value is 0.2? To see that LSQLIN is correct, you can generate contour plots of the function being minimized using the code below. You will see from the plots that as you descend through the isocontours, you converge precisely to the point on the boundary that LSQLIN is giving you.

    >> cc,cy

    cc =

      200.0000   27.5925
       27.5925   18.4254

    cy =

      107.8483
       11.7406

    [X,Y]=ndgrid(0.00001:.01:0.99998);
    xy=[X(:),Y(:)].';

    f=@(xy) sum((cc*xy).*xy)/2 -cy.'*xy; %objective function

    coef=[0.5392; 1e-5];  %alleged minimum
    fmin=f(coef);

    Z=reshape( f(xy)  , size(X));

    v=linspace(fmin+.05,fmin+5,10);
    contour(X,Y,Z,v); 
    title 'Contours of Objective Function'

could you also please let me know if I am using the fmincon code corresctly?

If you want to double check the results with FMINCON, you need to provide 2x1 initial guess. y will not be the correct size. A more sensible initial point would be C\y, so

 coef = fmincon(fun,C\y,[],[],[],[], lb ,ub,[],options);

[dav]

thank you very much Matt.

To double check, can I just use

coef = fmincon(fun,C\y,[],[],[],[], lb ,ub,[],options);

I think I have to write fun as a separate .m file.

I used the following two .m files but it gave me an error. am I doing it correct?

thanks

1)

function f = myfun(x)
f = norm(C*x(:)-y)^2;

2) (with the rst of the code)

            y = yt1;
            x= xt1;
            len = length(y);
            C = zeros(len,p);

            C(:,1) = 1; % The first column is for a0
            %for i = 1:p %Then we create shifted columns ( p in number ) for a 
             %  C(1+i:len,1+i) = x(1:len-i,1);     
            %end

            C = [C x];

             lb = [0.00001 0.00001]';
             ub = [0.99998 0.99998]';

             options =optimset('Display','off','LargeScale','off');

             coef = fmincon(@myfun,C\y,[],[]);

                  alpha0 = coef(1);
                  ar1 = coef(2);

[Matt J]

I think I have to write fun as a separate .m file.

No. You don't.

Making myfun a separate file should have worked fine though, as long as you don't really have all the code in myfun on 1 line.

You should develop a habit of sharing more complete information. Here, for instance, you haven't pasted the full error messages you are getting, making it very difficult to diagnose what's gone wrong. By not sharing what you're seeing, you're also being led to false conclusions, like the conclusion that you need myfun in a separate file. It's also been very hard to get the C.'*C etc... data from you. I'm just saying....

[Matt J]

In addition to my comments from 11 minutes ago, omiting any bounds here

 coef = fmincon(@myfun,C\y,[],[]);

will make comparison with LSQLIN (where you do apply bounds) invalid.

[dav]

I am really sorry for my mistakes. I did not see you had asked for C'C data for along time. I had to work on a different project and did not have time to ckeck on this.

I will be very careful in the future.

And I really appreciate all you valuable help.

Thank you very much!

 Here's what I used

              y = yt1;
              x= xt1;
              len = length(y);
              C = zeros(len,p);

              C(:,1) = 1; 
              C = [C x];

               lb = [0.00001 0.00001]';
               ub = [0.99998 0.99998]';

               options =optimset('Display','off','LargeScale','off');

               fun=@(x) norm(C*x(:)-y)^2;

               x = fmincon(fun,C\y,[],[],[],[], lb ,ub,[],options);

                    alpha0 = x(1);
                    ar1 = x(2);

But it gave me the following error

Error

Error: File: fmincon.m Line: 2 Column: 13
Unexpected MATLAB expression.

Error in aaa (line 147)
             x = fmincon(fun,C\y,[],[],[],[], lb ,ub,[],options);

Thanks.

[Matt J]

Error: File: fmincon.m Line: 2 Column: 13

That's weird, but in any case the error is coming from FMINCON itself, not from your objective function. Is it possible you edited fmincon.m by mistake? Copy/paste line 2 in fmincon.m for us so that we can see.

Also, I made a mistake when I said your separate file myfun is fine. You haven't passed C or y to myfun, so if you really want it in a separate file you would have to rewrite it

function f = myfun(x,C,y)
 f = norm(C*x(:)-y)^2;

and then do

 coef = fmincon(@(x) myfun(x,C,y),C\y,[],[]);

[dav]

This is what I see

fmincon('armarss',param0,A,b,Aeq,beq,LB,UB,'armanlc',options);

I did not write fmincon function. Do I have to write the fmincon function?

I thought I could use fmincon just the way I use lsqlin..

[Matt J]

No, the MathWorks wrote fmincon.m, but if you accidentally edited/corrupted it, you could have caused the kind of errors you see. I can't see anything wrong, but it might be a good idea to clean/reinstall MATLAB.

[dav]

I used the above code in another machine and it worked.

However, when I used in a different program (in my simulation program with a dbstop if naninf; code

it stopped at

coef = fmincon(fun,C1\y1,[],[],[],[], lb ,ub,[],options);

bug was in

function o = optimget(options,name,default,flag)

at

   value = defaultopt.(name);

the value was

value = inf

Do you have any idea what might be the problem?

thanks

Thi was what I got when I ran the program

NaN/Inf breakpoint hit for optimget>optimgetfast on line 108.
Stopping at next line.
108     value = defaultopt.(name);
109 end
K>> value

value =

     Inf

[Matt J]

Do you have any idea what might be the problem?

I doubt it's actually a problem. The MathWorks might set certain things to inf/nan in their code just to serve as default values.

[dav]

thanks.

I removed dbstop and re ran the code and the following warning was given

Warning: Rank deficient, rank = 2, tol =  5.130178e-14. 
> In fminseasonal at 710

at

bcoef = fmincon(fun,bC4\by4,[],[],[],[], lb ,ub,[],options);

Is it something to worry about?

Warning comes repeatedly in my simulation and tol = 5.130178e-14 changes

[Matt J]

Is fminseasonal your code? We need to know what matrix it is talking about when it says Rank deficient, rank = 2. Use "dbstop if warning" to stop there and see what matrix it is talking about.

[dav]

yes, fminseasonal is my code in which I do a simulation. When I use your command I got

Warning: Rank deficient, rank = 2, tol = 4.440892e-14. > In fminseasonal at 710

Warning from fminseasonal at 710
bcoef = fmincon(fun,bC4\by4,[],[],[],[], lb ,ub,[],options);

So I think it is talking about bC4 matrix. The difference between bC4 and the C4 matrix that you can see in my original code that I have posted here is that I use the coefficients I estimate to generate the data set again and do the same procedure.

Thanks a lot.

[Matt J]

Well, when this happens, check bC4 the same way we checked C. Do cond(bC4), show bC4.'*bC4, etc...

[dav]

ok thanks..

one last question,

in an earlier comment you had asked why I said that the correct value is 0.2.

I said it because the data set we used for that lsqlin code was generated such that the ar coeff was 0.2.

If lsqlin works properly shouldn't it give values closer to 0.2 rather than very small values which course all these issues?

Thank you very much for you time, answers and comments!!!!

[Matt J]

How close it will be to your simulated value depends on how much noise you have in your observation data.

Answer 2

[Youssef Khmou]

hi,

The determinant of the matrix bC is close to zero, you have to alter bC .

[dav]

how should i alter it please?

[Walter Roberson]

How did you create bC ? There might not be a way to avoid the problem for that set of data.

[dav]

I have a Time series data set (1000).

I took all the data corresponding to "day one" to one vector, day two to another vector and so on.. then I use lsqlin to estimate the parameters. However, the parameter estimates sometime are VERY CLOSE TO ZERO, but the actual parameters aren't that close to zero. this is why I am getting this warning.

Is there any way to correct it please?

CODE:

clc;
clear;
%warning off;

runs =100;

xxtt=[];
bxxtt=[];

p=1;

       yt1=[];
       yt2=[];
       yt3=[];
       yt4=[];
       yt5=[];

       xt1=[];
       xt2=[];
       xt3=[];
       xt4=[];
       xt5=[];

epsi=zeros(3000,1);
simsig=zeros(3000,1);

      a0(1)=0.5; a1(1)=0.2;

      for i=2:4
      a0(i)=0.1; a1(i)=0.3;
      end

      a0(5)=0.1;a1(5)=0.2;

for i = 1:3000

      m = mod(i,5);
      if  m==0
          m=5;
      end
                   if (i==1) 
                       simsig(i) = a0(m)/(1-a1(m));
                       s=(simsig(i))^0.5;
                       epsi(i) = normrnd(0,1) * s;

                   else                     
                       simsig(i) = a0(m)+ a1(m)*(epsi(i-1))^2;
                       s=(simsig(i))^0.5;
                       epsi(i) = normrnd(0,1) * s;
                   end
  end

epsi2=epsi.^2;

yt = epsi2(2001:3000);

ytinitial = epsi2(2001:3000);
sig = simsig(2001:3000);

ytlast=epsi(3000)^2;
ytlast2=epsi(2999)^2;

s=1;

for i=1:1000

       m = mod(i,5);
      if  m==0
          m=5;
      end

      if  m==1
          yt1(s) = yt(i);
          if i==1
          xt1(s)=0.00001;
          else
          xt1(s) = yt(i-1);
          end
      elseif m ==2
          yt2(s) = yt(i);
          xt2(s) = yt(i-1);
      elseif m==3
          yt3(s)= yt(i);
          xt3(s) = yt(i-1);
      elseif m==4
          yt4(s)= yt(i);
          xt4 (s)= yt(i-1);
      elseif m==5
          yt5(s)= yt(i);
          xt5(s)= yt(i-1);
      end

      s=s+1;

end

yt1(yt1==0) = [];
yt2(yt2==0) = [];
yt3(yt3==0) = [];
yt4(yt4==0) = [];
yt5(yt5==0) = [];

xt1(xt1==0) = [];
xt2(xt2==0) = [];
xt3(xt3==0) = [];
xt4(xt4==0) = [];
xt5(xt5==0) = [];

yt1=yt1';
yt2=yt2';
yt3=yt3';
yt4=yt4';
yt5=yt5';

xt1 = xt1';
xt2 = xt2';
xt3 = xt3';
xt4 = xt4';
xt5 = xt5';

ytt = [ yt1, yt2, yt3, yt4, yt5];
xtt = [ xt1, xt2, xt3, xt4, xt5];

%xxtt=[xxtt;xtt];

for m = 1:5

      y = ytt(:,m);
      x= xtt(:,m);
      len = length(y);
      C = zeros(len,p);

      C(:,1) = 1; % The first column is for a0
      %for i = 1:p %Then we create shifted columns ( p in number ) for a 
       %  C(1+i:len,1+i) = x(1:len-i,1);     
      %end

      C = [C x];

       lb = [0.00001 0.00001]';
       ub = [0.99998 0.99998]';

       options =optimset('Display','off','LargeScale','off');

       coef = lsqlin(C,y,[],[],[],[], lb ,ub,[],options);

      alpha0(m) = coef(1);
      ar1(m) = coef(2);
  end