Improve performance of code to solve equation

Question

I am quite novice to MATLAB and I am struggling to find a solution to this equation.

Where the matrices' dimensions are {λ} = N×K, {Y} = N×D, {π} = 1×K, and {μ} = D×K.

What I have created looks more like a monstrosity than an efficient MATLAB code, and that is due to a serious lack of MATLAB skills.

Actually, I decided to solve this step by step, in order to get myself familiarised with how MATLAB works. I ended up with an ultra inefficient code that I cannot even combine it. Even the final result is wrong, as I want to create a NxK matrix. Any guidance would be much appreciated.

%Dimensions:
nn = 10;
dd = 7;
kk = 5;

%Initial variables:
lambda0 = rand(nn,kk);
sigma = rand(1);
Y = rand(nn,dd);
mu = rand(dd,kk);

%Calculate pies:
First_part = log(pie(:)./(1.-pie(:))); % Kx1 vector

%Calculate Second part:
for n = 1:nn
for d = 1:dd
lambda_mu(d,:) = lambda0(n,:).*mu(d,:); %lambdamu is a DxK matrix
end
lambda_mu2(n,:) = sum(lambda_mu,2); %This is a NXD matrix
end

%Y-lambdamu2:
for n = 1:nn
YY(n,:) = Y(n,:)-lambda_mu2(n,:); % This is a NxD vector
end

%YY*mu:
Second = (YY*mu)./sigma^2; % NxK

%Third:
Third = (mu'*mu)./2*sigma^2; % KxK

%Final:
for n=1:nn
Final_part = transpose(First_part(:))+Second(n,:)-Third;
end

Update1: Ok, I owe myself a beer. So, I made some progress to create a {λ} = N×K matrix by modifying my code as follows:

Update[2]: Spending too much time programming and you miss details. I fixed the error as I had set {Y} = D×K, instead of {Y} = D×N.

%Dimensions:
nn = 10;
dd = 7;
kk = 5;

%Initial variables:
lambda0 = rand(nn,kk);
sigma = rand(1);
Y = rand(dd,nn);
mu = rand(dd,kk);
pie = rand(1,kk);
sigma = rand(1);

%Calculate the equation.
for n = 1:nn
for k = 1:kk
    lambda(n,k) = log(pie(k)/(1-pie(k))) + (transpose(Y(:,n)-... 
        (sum(lambda0(n,1:end ~= k).*mu(:,1:end ~= k),2)))*mu(:,k))/...
        sigma^2 - transpose(mu(:,k))*mu(:,k)/2*sigma^2;
end
end

But, now the problem is that it goes only up to n=5, when I try to loop for n=6 for example I get this message: "Index exceeds matrix dimensions", so practically I only get a 5×5 matrix. Suggestions please?

P.S: I even tried to change the lambda from lambda(n,k) to lambda(k,n) but yet the result is the same.

Answer 1

In order to get dimensions working, multiply by ones(nn,kk) and ones(1,kk)

nn = 10;
dd = 7;
kk = 5;

%Initial variables:
lambda0 = rand(nn,kk);
sigma = rand(1);
Y = rand(nn,dd);
mu = rand(dd,kk);
pie = rand(kk,1);

%Calculate pies:
First_part = ones(nn,kk)*log(pie./(1.-pie))*ones(1,kk); % Kx1 vector

%Calculate Second part:
for n = 1:nn
    for d = 1:dd
        lambda_mu(d,:) = lambda0(n,:).*mu(d,:); %lambdamu is a DxK matrix
    end
    lambda_mu2(n,:) = sum(lambda_mu,2); %This is a NXD matrix
end

%Y-lambdamu2:
for n = 1:nn
    YY(n,:) = Y(n,:)-lambda_mu2(n,:); % This is a NxD vector
end

%YY*mu:
Second_part = (YY*mu)./sigma^2; % NxK

%Third:
Third_part = ones(nn,kk)*(mu'*mu)./2*sigma^2; % KxK


%Final:
for n=1:nn
    Final_part(n,:) = First_part(n,:)+Second_part(n,:)-Third_part(n,:);
end

EDIT: I was thinking more like:

%Initial variables:
lambda = rand(nn,kk);
sigma = rand(1);
y = rand(nn,1);
mu = rand(nn,kk);
pie = rand(1,kk);

%Calculate pies:
First_part = ones(nn,1)*log(pie./(1-pie)); % NxK vector

to avoid loops.

Maybe use something like:

ones(kk)-diag([ones(kk,1)])

when doing the sum over i not equal to j, or rather use a third dimension to represent j.