TLDR: I am trying to do maximum likelihood fitting of a dataset having two mixed populations, observed over a subset of their parameter space, within it to two pdfs. I include working code with gaussian examples.
Hello, I am new to this site. My question (and example code) is also long, so I will host some information in pastebin. I apologize, I am neither a mathematician nor a computer scientist, so some things may be ugly or inefficient.
I have a set of data which has samples of two populations within it (we have a good prior of what percentage of the data belongs to each population). Neither sample probes the full spatial extent of the data, but the probability density function (pdf) is presumed to be smooth. I use a wide variety of pdfs, but for this discussion I will utilize only Gaussians for simplicity-- however, I manually normalize the Gaussians, since my pdfs are not normalized.
We can think of the data as being sampled from two different Gaussians, with some kind of cut off boundaries, beyond which the population continues, but the sampling ends. The boundaries of the two populations intersect for a range.
The metacode would be something along the lines of:
A) Make a dataset by sampling two pdfs with known parameters, within two different given ranges.
$$ pdf = gaussian = e^{\frac{-1}{2}(\frac{x-mu}{sigma})^{2}} $$
B) Guess some parameters
--1) Find the normalizing constants ($n_{1}$ and $n_{2}$) for these parameters by integrating the Gaussians of this guess analytically over the observed space. $$ n = \sum{volume*density} $$ $$ density = pdf(x) :: volume = step size $$
--2) For each item in the dataset: $$ Probability (p)= \frac{ratio}{n_{1}}*pdf_{1}(item) + \frac{1-ratio}{n_{2}}*pdf_{2}(item) $$
--3) Sum the log likelihoods for all data to find the total likelihood of the guess $$ L = \sum{log(p)} $$
C) Use some likelihood landscape explorer to converge on a best answer (here I use brute grid for simplicity)
I have had great success disentangling two populations if they are sampled in their entirety, but poor success disentangling two populations that have been truncated (the answers are frequently nonsensical and boundary dominated).
The code at (~100 lines, but not dense code) will run as is. It follows the metacode above and successfully finds the parameters of the two Gaussians sampled over (essentially) their entire area. However, changing lines 31 and 32 to:
range1=[-2,6]
range2=[2,12]
for example, causes the algorithm to find the incorrect values. My question is why.
Edit: Thank you for the suggestion, I've reformatted this to be easier to read.