0

I was wondering if it is possible to use NMinimize from mathematica with an objective function, which contains random variables? E.g. I have a function with parameters which follow a distribution (normal and truncated normal). I want to fit its histogram to data that I have and constructed an objective function which now I need to minimize (so, the objective function depends on the mus and sigmas of the parameters and need to be determined). If I run my code, there is an error message: It claims the parameter for the NormalDistribution needs to be positive (If I plug in numbers for the mus and sigmas of my objective functionby hand, i don't get an error message). So, I am wondering if it is not possible for NMinimize to handle a non-analytic function. Thanks! Here, I give you an example code (please note that the original function is more complicated)


listS and listT are both lists of event times. I want to fit the curve of my statistical model for the times (here, a very simple one, it consists of a truncated normal distribution) to the data I have. For this I compare the survival curves and need to minimize the sum of the least squares. My problem is that the function NMinimize doesn't seem to work. (Please note, that the original objective function consists of a more complicated function with parameters that are random variables)

(* Both lists are supposed to be the list of times *)

SurvivalS[listeS_, x_] := Module[{res, survivald},
survivald = SurvivalDistribution[listeS];
res = SurvivalFunction[survivald, x];
res]
Residuum[listeT_, listeS_] := 
Table[(SurvivalS[listeT, listeT[[i]]] - SurvivalS[listeS, listeT[[i]]]), {i,
 1, dataN}];
LeastSquare[listeT_, listeS_] := 
Total[Function[x, x^2] /@ 
Residuum[listeT, 
 listeS]];(* objective function, here ist is the sum of least square *)

objectiveF[mu_, sigma_] := 
Piecewise[{{LeastSquare[listeT, listeS[mu, sigma]], mu > 0 && sigma > 0}}, 
20 (1 + (sigma + mu)^2)];

pool = 100; (* No. points from MonteCarlo *)

listeS[mu_, sigma_] := RandomVariate[TruncatedDistribution[{0, 1}, NormalDistribution[mu, sigma]],pool];(* simulated data *)
listeT = Sort[RandomVariate[TruncatedDistribution[{0, 1}, NormalDistribution[.5, .9]],60]]; (* list of "measured" data *)
dataN = Length[listeT];
NMinimize[objectiveF[mu, .9], {{mu, .4}}]

The error message is: "RandomVariate::realprm: Parameter mu at position 1 in NormalDistribution[mu,0.9] is expected to be real. >>"

  • There is a dedicated site for [mathematica.se], and most of us have moved over there. If you'd like, we can move the question over, too. – rcollyer Aug 25 '12 at 01:16
  • NMinimize[objectiveF[mu, .9], {{mu, .4}}] is not the right structure. Try something like NMinimize[objectiveF[mu, .9], {mu, .2,1.2}]. Also it might be necessary to restrict objectiveF to only numeric arguments, via objectiveF[mu_?NumberQ...] := ... Finally you might look into memoizing objectiveF since it might be relatively expensive to evaluate. Could also use Compile in one or more places for further speed boost. – Daniel Lichtblau Aug 26 '12 at 21:32
  • rcollyer: yes please, move it to the dedicated mathematica site! daniel: thank you for your comment, i will try that... – user1620246 Sep 02 '12 at 12:01

0 Answers0