From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: Alternative to matrix inverse?
Date: Thu, 23 Feb 2012 23:47:20 -0800 (PST)
Lines: 105
Message-ID: <>
References: <jhvcv1$j70$> <ji0bpa$m1k$>
 <ji0ehf$28d$> <> <ji722q$d29$>
Mime-Version: 1.0
Content-Type: text/plain; charset=windows-1252
Content-Transfer-Encoding: quoted-printable
X-Trace: 1330071117 3007 (24 Feb 2012 08:11:57 GMT)
NNTP-Posting-Date: Fri, 24 Feb 2012 08:11:57 +0000 (UTC)
Injection-Info:; posting-host=; posting-account=X3eThQoAAACh5vOSip_rNRzKAq7k0jPW
User-Agent: G2/1.0
X-Google-Web-Client: true
X-Google-Header-Order: ARLUEHNKC
X-HTTP-UserAgent: Mozilla/4.0 (compatible; MSIE 8.0; Windows NT 5.1;
 Trident/4.0; .NET CLR 2.0.50727; .NET CLR 3.0.04506.30; .NET CLR
 3.0.04506.648; .NET CLR 3.0.4506.2152; .NET CLR 3.5.30729; .NET4.0C; .NET4.0E),gzip(gfe)
Xref: comp.soft-sys.matlab:758818

On 24 Feb., 05:01, "Jean-Philippe " <> wrote:
> Christopher Creutzig <> wrote in message <>...
> > On 21.02.12 16:51, Jean-Philippe wrote:
> > > Thank you very much for your replies.
> > > In fact, P is a graph representing the weights of transitions between the neighbouring pixels for every pixels. The inverse is then used to extend que probabilities to all the other pixels. We obtain B, which contains probabilities for every pixels. The problem with the iterative approach I suggested is that the "spreading" of the probabilities accross B would require too many iterations and would likely be slower than inv.
> > I.e., you are using
> >   infinity
> >      ---
> >      \              i        1
> >      /     (alpha P)  = -----------
> >      ---                1 - alpha P
> >     i = 0
> > to find the sum on the right hand side, as the limit of some Markov
> > chain, right? (If the formula doesn't display correctly, try a
> > monospaced font. That may seem outdated, but continues being the best
> > choice for Usenet, really …)
> > > Is there another way to propagate the probabilities accross B without using the inverse, or is there a way to parallelize this equation?
> > You could simply calculate the sum directly, truncating it at some
> > modest upper limit, say, by monitoring norm((alpha*P)^i). If
> > norm(alpha*P) is sufficiently small, the sum should converge rather
> > quickly – plus, you avoid potential cancellation errors on the diagonal
> > for really small entries.
> > Christopher
> Hi Christopher, thank you for helping me out :)
> The problem I am actually trying to solve is to use the matrix B as the probabilities of transition between every two pixels of the image. To find B, we use the weight matrix P which contains the intensity difference for every pixels and their neighbouring pixels. We use alpha to control how the probabilities are dispersed across B. Those probabilities are then mixed with others seed pixels probabilities and we achieve image segmentation by keeping the most probable seed class for every pixel.
> I am still learning all those maths right now, so I am not sure what is a Markov Chain *blush*...
> What I was asking first, was a way to obtain a similar result than what the function inv gives me because it is really slow when I'm working on big matrices. P and B are NxN matrices where N is the total number of pixels.
> I suggested an iterative method to estimate B, but it requires many iterations to obtain something useful and it becomes even slower than inv. I was wondering if I could optimize the matlab implementation of this method (M=M*P; C=C+M;) where all the matrices are currently sparses. Therefore I cannot use the function norm on those matrices and I cannot do full() on any of those because it causes an OUT OF MEMORY error.
> I tried my best to understand Roger's and your suggestions, but it seems that I don't get it... Thank you for your patience with my problem.
> What I understood from Roger's anwser is that I could modify the iterative method to be more efficient by reducing the number of operations. I have been able to implement a way to extract all the diagonals of matrix P in a 9xN matrix. But then, I do not know what to do next to iterate.
> From your answer, I understand that I could use a third way to calculate B by iterating (alpha*P)^i and by looking at the norm of P to stop the iterations automatically based of a certain value. If I understood correctly, then how can I get the norm of a sparse matrix?
> --- as reminder ---
> using the matlab inv
>     B = ((1-alpha) * inv(speye(N)-(alpha*P)));
> or the iterative way
>     B = (1-alpha) * inverse(speye(N), (alpha*P), limit);
>     B = (1-alpha) * C;
>     C = I;
>     for i = 1:limit
>         I = I*P;
>         C = C+I;
>     end- Zitierten Text ausblenden -
> - Zitierten Text anzeigen -

A possible norm would be the maximum norm, i.e. the maximum of the
absolute values of the matrix
But for the infinite series to converge you must have that the
spectral radius of the matrix alpha*P
is less than 1. Is this the case for your application ? If you don't
know, you could make a test and check
whether (alpha*P)^i is near to the zero matrix for large values of i.

Best wishes