Re: Speedup/Refaktoring loop? Extended...

[Jaroslav Hajek]

On Tue, Aug 19, 2008 at 2:36 PM, Andreas Romeyke <address@hidden> wrote:
> Hello,
>
> Am Tue, 19 Aug 2008 13:20:42 +0200
> schrieb "Jaroslav Hajek" <address@hidden>:
>
>>
>> try:
>>
>> for iteration = 1:maxit
>>    a=transpose(U) * X;
>>    tmpU = transpose(U)*U;
>>    for i=1:n
>>        bi=tmpU(i,:)*V;
>>        V(i,:) = V(i,:) .* a(i,:) / b;
>>    end
>>
>>    .. and analogously for the second loop
>> end
>
> Thanks to Jaroslav the algorithm could be described much simpler.
> Because the paper of Lee and Seung is not very clear the update rules
> seems more difficult as exist.
>
> In the paper they wrote:
>
>              (U^TX)_ij
> V_ij <-- V_ij ---------------
>              (U^TUV)_ij
>
> and
>              (XV^T)_ki
> U_ki <-- U_ki ---------------
>              (UVV^T)_ki
>
> I thought that every rule infects the other one by incrementing indices
> k,i and j. But in reality (tested in experiments) the stuff could be
> written in terms of Octave as:
>
>    a=transpose(U) * X;
>    b=transpose(U) * U * V;
>    V=V .* a ./ b;
>    a=X * transpose(V);
>    b=U * V * transpose(V);
>    U=U .* a ./ b;
>

Well, that's what we guessed from the beginning. If it was the other
way around, it would be quite an inefficient update. But the paper
really is unclear; it requires some guessing to get it right.
A final remark: why do you use transpose() instead of the ' and .'
operators? Octave 3.2.x will optimize expressions like U.'*U or X*V.',
but not if you use transpose().

> Bye Andreas
> --
> Software Developer / Dipl. Inform. (FH)
> Max Planck Institute for Human Cognitive and Brain Sciences
> Department of Psychology
> Stephanstr. 1a, 04103 Leipzig, Germany
>
>



-- 
RNDr. Jaroslav Hajek
computing expert
Aeronautical Research and Test Institute (VZLU)
Prague, Czech Republic
url: www.highegg.matfyz.cz