# Fitting Gaussian to data: Density-Estimation vs Regression

Recently I discussed the following topic with a friend. The setting is that we have a one-dimensional set of data. (In the example it was points of students, we would be grading.) The goal was to make a density estimation, but not to use anything "fancy" like kernel density estimation, but just use a Gaussian as the estimation. (Sure that makes the big assumption that the data is Gaussian, but that is not the point here.) We discussed two ways:

1. Make a density estimation using an unsupervised learning method, e.g. using EM-algorithm. In this case the claim is that simply calculating the mean of the data and the standard derivation is already giving one the parameters to get the Gaussian parametrized the right way.
2. Add up the number of occurrences for each value and then use a supervised learning regression with the Gaussian as function, powered by an optimisation algorithm.

Through discussion we found out that the two clearly have different outcomes - though probably similar. For case 2 we optimise the parameters of the Gaussian such, that the sum of the distances from the Gaussian to the occurrences is minimized (along the y-axis if you will). For case 1 we optimise the parameters along the x-axis if you will.

## Questions

1. Preliminary question: Is it correct that the EM-algorithm will have the same result as just calculating mean and standard deviation from the data?

Assuming the answer to the preliminary question is "yes, the results are the same":

1. What is the intuition interpretation of those two approaches?
2. None of them can be wrong in itself, but one of them could be wrong in the sense of wrong-usage. Meaning: I want to do something and I use the wrong method for it, because of misunderstood interpretation of what happened. So in that sense: Is one of them wrong in a way?

## Example Code

I managed to express myself in R code. As one can see from the plots, the result is definitely not the same for any dataset. Only if the data is Gaussian and large the results get to be similar. But that means little to me ($2^x$ and $exp(x)$ converge both to $infinity$ for $x -> infinity$ and they have little in common otherwise).

library(ggplot2)
library(dplyr)
library(tidyr)

std <- 10

datasets <- list(
data.frame( x = round(rnorm(1000, sd = std))),
data.frame( x = round(rnorm(1000, sd = std))) %>% filter(x > 0)
)

for(data1 in datasets){

data1$densest <- dnorm(data1$x, mean = mean(data1$x), sd = std, log = FALSE) data2 <- data1 %>% group_by(x) %>% summarise(count = n()) f <- function(x, m, sd, k) { k * exp(-0.5 * ((x - m)/sd)^2) # 1/sqrt(2*pi*sd^2) * } cost <- function(par) { rhat <- f(data2$x, par, par, par)
sum((data2$count - rhat)^2) } o <- optim(c(0, std, 10), cost, method="BFGS", control=list(reltol=1e-9)) data1$regr <- f(data1$x, o$par, o$par, o$par)
data1$regrNormalized <- dnorm(data1$x, mean = o$par, sd = abs(o$par), log = FALSE)

g1 <- ggplot(data=data1, aes(x=x)) +
geom_point(data=data2, aes(x=x, y=count), alpha=0.2)+
geom_line(aes(y=densest), color="green") +
geom_point(data = data.frame(mean = mean(data1$x)), aes(x=mean, y=0), color="green") + geom_line(aes(y=regr), color="blue") + geom_point(data = data.frame(mean = o$par), aes(x=mean, y=0), color="blue")

g2 <- ggplot(data=data1, aes(x=x)) +
geom_line(aes(y=densest), color="green") +
geom_point(data = data.frame(mean = mean(data1$x)), aes(x=mean, y=0), color="green") + geom_line(aes(y=regrNormalized), color="blue") + geom_point(data = data.frame(mean = o$par), aes(x=mean, y=0), color="blue")

plot(g1)
plot(g2)

}


### Plots

• What's wrong with the maximum-likelihood estimators of the mean and variance of your Gaussian given the data? ie sample mean and mean square deviation from the sample mean? Closed form solutions, no messing. Or are you really interested in comparing outcomes of different methods for different (less tractable) distributions? Dec 2, 2016 at 9:20
• @Spacedman: Yes I am. It is an actual data science question, not a practical, applied, works-so-I-am-fine-with-it business analyst question. (Apologies, for my cynicism - too many data non-scientists around theses days. You know what I mean ;-).) This is about the intuitive interpretation of the mathematically different approaches. It's a math question. Dec 2, 2016 at 18:49
• If I understand method (2) correctly, you are binning the data first into counts - this loses information. That should be enough for intuition to see that you could get different results - or at least you are estimating different things. Dec 4, 2016 at 14:26
• @Spacedman: No, it does not lose information. The set {1,1,1,1,2,2,2,3,4,4,4,4,4,5,5} equals the map {1->4, 2->3, 3->1, 4->5, 5->2} in information. The reason those are different because of different loss functions (where the loss function in 1 is implicit, just as the loss function for K-Means is also just implicit - as an example.) Dec 5, 2016 at 15:16
• ah, I was thrown by your Gaussian assumption again. A set of integers can't be from a Gaussian because continuous. But you don't care about that. If you could include some code for what you did it might be possible to better see what you are doing and if the asymptotic results (for large amounts of data) converge - because its still a bit unclear. Dec 5, 2016 at 15:24