I need to fit a non-analytic data set to another data set using at least 2 parameters: shift and scale. For my 2 data sets (data1 and data2), I will to apply these 2 parameters to data1 such that it matches data2. This works perfectly for only one of those parameters, but when I try to do both at the same time, it fails. The 2nd parameter is never actually fit...
I will illustrate this by creating a "double hump" distribution using two normal distributions.
Clear["Global`*"]
NumPerSet = 1000;
(* This is the first "double-hump" distribution *)
set1 = RandomVariate[NormalDistribution[6, 1.2], NumPerSet];
set2 = RandomVariate[NormalDistribution[11, 1.2], NumPerSet];
data1 = Join[set1, set2];
Then I do the same thing again, but with the peaks shifted to the left by 1:
(* This is like data1 but "shifted" to the left by 1 *)
set3 = RandomVariate[NormalDistribution[5, 1.2], NumPerSet];
set4 = RandomVariate[NormalDistribution[10, 1.2], NumPerSet];
data2 = Join[set3, set4];
Then, I apply the scale factor, where each data point in the set has been reduced by 10%:
(* Now apply a scale factor to the data2 set *)
data2 = data2 0.9;
To visualize this, I can plot:
(* Create some plottable distributions *)
skd1 = KernelMixtureDistribution[data1];
skd2 = KernelMixtureDistribution[data2];
(* Show each distribution and the associated histogram *)
Show[Plot[{PDF[skd1, x], PDF[skd2, x]}, {x, 0, 15}, Frame -> True, PlotLabel -> "The Setup"], Histogram[{data1, data2}, Automatic, "PDF", ChartLegends -> {"data1", "data2"}]]
Which returns:
Great, now to solve for the parameters:
(*Mean of distribution*)
\[Mu] = Expectation[x, x \[Distributed] skd2];
(*Estimate a by method of moments*)
mom = Abs[Mean[data1] - \[Mu]];
(*Estimate a by maximum likelihood*)
logL = LogLikelihood[skd2, (data1 - shift) scale];
fmout = FindMaximum[logL, {{shift, mom}, {scale, .9}}, AccuracyGoal -> 20]
Print["fmout = ", fmout]
Print["Max Likelihood value for:"]
Print[" Shift = ", Part[Part[Part[fmout, 2], 1], 2]]
Print[" Scale = ", Part[Part[Part[fmout, 2], 2], 2]]
Print[" MLE = ", Part[fmout, 1]]
(* Now show the actual adjustment on data1 which supposedly fits data2 the best *)
Show[Plot[{PDF[skd2, x]}, {x, 0, 15}, Frame -> True, PlotLabel -> "Fitted Values"], Histogram[{data2, (data1 Part[Part[Part[fmout, 2], 2], 2]) - Part[Part[Part[fmout, 2], 1], 2]}, Automatic, "PDF", ChartLegends -> {"data1", "data2"}]]
(* Now show what we already know is the true shift and scale *)
Show[Plot[{PDF[skd2, x]}, {x, 0, 15}, Frame -> True, PlotLabel -> "Already Known Values"], Histogram[{data2, (data1 - 1) .9}, Automatic, "PDF"]]
What consistently comes back is that the first parameter (shift) is fit, but the second parameter is just stuck at whatever I put in for the starting value (in this case, 0.9).
For example, in one case, the shift was 1.76895, scale was 0.9, and the plot:
And by applying the "already known" values of shift = 1, scale = 0.9, I get:
Obviously a much better fit... Can anyone help me understand why this is failing?
Also, if I try to plot the MLE via ContourPlot, I get something that looks like the whole thing is just borked (shift on X-axis, scale on Y-axis, red dot indicates the "fitted" parameters):
Thanks!
edit: Code for the contour plot:
Show[ContourPlot[logL, {shift, 0, 2}, {scale, 0.8, 1.0}], ListPlot[{{Part[Part[Part[fmout, 2], 1], 2], Part[Part[Part[fmout, 2], 2], 2]}}, PlotStyle -> Red]]
NMaximize[{logL, 0.5 <= shift <= 1.5, 0.8 <= scale <= 1.0}, {shift, scale}]which yielded similarly bad results. This time it just sets shift as low as possible... $\endgroup$LogLikelihoodsince it doesn't seem to work with multiple parameters (for non-pre-determined distributions), but this gives decent results:nmaxOut = NMaximize[{DistributionFitTest[(data1 - shift) scale, skd2, "CramerVonMises"], 0 <= shift <= 1.5, 0.8 <= scale <= 1.0}, {shift, scale}]$\endgroup$