Which is the best optimization method for optimizing over probability mass functions?

I have a nonlinear scalar function f of a n-dimensional variable/vector p, where p is a probability mass function (PMF) for a (discrete) random variable with finite support. This means that the entries of the vector p are probabilities that sum to one, so and for .

What's the best MATLAB optimization method for finding the minimum of f? I have been using fmincon,but I thought perhaps there's a better method.

Here's my code for a random toy nonliner function f (which is the function funOpt in the code):

%Optimization problem  for a probability mass function (PMF)
numbPMF=3; %number  probability mass function (PMF) values
rng(1); %set random seed for reproducibility
p0=rand(1,numbPMF);p0=p0/sum(p0); %create a random PMF

funOpt=@(p)(sum((p-p0).^2+(p-p0).^3)); %example of a nonlinear function
pmf0=zeros(1,numbPMF); %initial guess of PMF

probTotal=1; %total probability (default is one)

%%%START Optimization parameters START%%%
%For details, see https://au.mathworks.com/help/optim/ug/fmincon.html
lb = zeros(1,numbPMF); %lower bounds for probabilities
ub = ones(1,numbPMF); %upper bounds for probabilities
nonlcon = @(x)funSumProb(x,probTotal); %nonlinear constraint; see below

%leave all the other optimization parameters void
A = []; %linear equality contraint ie A*x = b
b = []; %linear contraint ie A*x =b
Aeq = []; %linear inequality contraint ie A*x <= b
beq = []; %linear inequality contraint ie A*x <= b
%%%END Optimization parameters END%%%

[pmfMin,F]=fmincon(funOpt,pmf0,A,b,Aeq,beq,lb,ub,nonlcon);

p0
pmfMin
errorPMF=mean(abs(pmfMin-p0))

%define function for nonlinear constraint
%contraint is the probabilities summing to one
function [c,ceq] = funSumProb(x,totalSum)
ceq= (sum(x)-totalSum); %sums to one
c = [];
end

Thanks in advance.

Answer 1

[Matt J]

fmincon is most appropriate because your objective is non-linear and non-quadratic. However, your constraints are all linear, so you shouldn't be using nonlcon:

Aeq = ones(1,numbPMF); 
beq=1;

[pmfMin,F]=fmincon(funOpt,pmf0,[],[],Aeq,beq,lb,ub,[]);

[H. Paul Keeler]

Right. That makes sense. Thanks very much for that.