Separating two populations from the sample I'm trying to separate two groups of values from a single data set. I can assume that one of the populations is normally distributed and is at least half the size of the sample. The values of the second one are both lower or higher than the values from the first one (distribution is unknown). What I'm trying to do is to find the upper and lower limits that would enclose the normally-distributed population from the other.
My assumption provide me with starting point:


*

*all points within the interquartile range of the sample are from the normally-distributed population.


I'm trying to test for outliers taking them from the rest of the sample until they don't fit into the 3 st.dev of the normally-distributed population. Which is not ideal, but seem to produce reasonable enough result.
Is my assumption statistically sound? What would be a better way to go about this?
p.s. please fix the tags someone.
 A: *

*For data in IQR range you should use
truncated normal distribution (for
example R package gamlss.tr) to
estimate parameters of this
distribution.  

*Another approach is using mixture models with 2 or 3 components (distributions). You can fit   such models using gamlss.mx package (distributions from package gamlss.dist can be specified for
each component of mixture).

A: This assumes that you don't even know if the second distribution is normal or not; I basically handle this uncertainty by focusing only on the normal distribution. This may or may not be the best approach.
If you can assume that the two populations are completely separated (i.e., all values from distribution A are less than all values from distribution B), then one approach is to use the optimize() function in R to search for the break-point that yields estimates of the mean and sd of the normal distribution that make the data most likely:
#generate completely separated data
a = rnorm(100)
b = rnorm(100,10)
while(!all(a<b)){
    a = rnorm(100)
    b = rnorm(100,10)
}

#create a mix
mix = c(a,b)

#"forget" the original distributions
rm(a)
rm(b)

#try to find the break point between the distributions
break_point = optimize(
    f = function(x){
        data_from_a = mix[mix<x]
        likelihood = dnorm(data_from_a,mean(data_from_a),sd(data_from_a))
        SLL = sum(log(likelihood))
        return(SLL)
    }
    , interval = c(sort(mix)[2],max(mix))
    , maximum = TRUE
)$maximum

#label the data
labelled_mix = data.frame(
    x = mix
    , source = ifelse(mix<break_point,'A','B')
)
print(labelled_mix)

If you can't assume complete separation, then I think you'll have to assume some distribution for the second distribution and then use mixture modelling. Note that mixture modelling won't actually label the individual data points, but will give you the mixture proportion and estimates of the parameters of each distribution (eg. mean, sd, etc.).
A: If I understand correctly, then you can just fit a mixture of two Normals to the data. There are lots of R packages that are available to do this. This example uses the mixtools package:
#Taken from the documentation
library(mixtools)
data(faithful)
attach(faithful)

#Fit two Normals
wait1 = normalmixEM(waiting, lambda = 0.5)
plot(wait1, density=TRUE, loglik=FALSE)

This gives:
Mixture of two Normals http://img294.imageshack.us/img294/4213/kernal.jpg
The package also contains more sophisticated methods - check the documentation.
A: I'm surprised nobody suggested the obvious solution:
 #generate completely separated data
library(robustbase)
set.seed(123)  
x<-rnorm(200)
x[1:40]<-x[1:40]+10  
x[41:80]<-x[41:80]-10
Rob<-ltsReg(x~1,nsamp="best")
#all the good guys
which(Rob$raw.weights==1)

Now for the explanation: the ltsReg function 
in package robustbase, when called with the option 
nsamp="best"

yields the univariate (exact) MCD weights. 
(these are a n-vector 0-1 weights stored 
in the $raw.weights object. The algorithm
 to identify them is the MCD estimator (1)).  
In a nutshell, these weights are 1 for the 
members of the subset of $h=\lceil(n+2)/2\rceil$ most 
concentrated observations. 
In dimension one, it starts by sorting all the 
observations then computes the measure of all 
contiguous subsets of $h$ observations: denoting 
$x_{(i)}$ the $i^{th}$ entry of the vector of sorted 
observations, it computes the measure of
(e.g. $(x_{(1)},...,x_{(h+1)})$ then $(x_{(2)},...,x_{(h+2)})$ 
and so forth...) then retains the one with smaller 
measure.
This algorithm assumes that your group of interest numbers
 a strict majority of the original sample and that it
has a symmetrical distribution (but there no
 no hypothesis on the distribution of the remaining 
$n-h$ observation).

(1) P.J. Rousseeuw (1984). Least median of squares regression, Journal of
  the American Statistical Association.

