In reply to Rick Wicklin.

Rick,
First of all thanks for your time and for the answer. I have to ask more to understand better. Here is what I did.

I made a small test to check the CI with MCSEM (as you said), I run the simulation for 1000 times then count means between CI for 95%. 948 of them lied in the interval. When I try for the second interval I gave above (as you interpret as wrong) 53 of the values lied in the interval.

So, roughly I can think this as we have a population (unknown), we take 1000 samples with sample size 10.000 for each. Right?

I run the simulation for one time (as you said) with 50.000 Dots. So, I have 1 sample with size 50.000. Let's call this case 1:

Mean 3.1287
Standard Error 0.0001
Median 3.1355
Mode 3.1111
Standard Deviation 0.0226

Then I run it 1000 times for 50.000 Dots. This is case 2:
Mean 3.1417
Standard Error 0.0002
Median 3.1418
Mode 3.1423
Standard Deviation 0.0071

Here is my question. If I understand correctly, SEM in case 1 (0.0001) is an estimate of the SD in case 2 (0.0071). This should lead me to CI(case 2) = (3.1278, 3.1556) (Mean +/- 1.96 * SD) which already covers 94.8% of the data.
And for case 1, by using SEM CI(case 1) = (3.1274, 3.1290) (Mean +/- 1.96 * SEM * 4(?**)) which only covers 0.19% of the data in case 2. What I am missing here? Obviously, 0.19% is not a good approximation, I am doing something wrong.

?** 4 didn't change much here, since SEM is 0.0001. Values I typed are without 4.

By the way, the distribution of case 1 (1 sample 50.000 dots) is not normal. (Kurtosis 267.456, Skewness 7.429) Is this related? Do I have to go for case 2 since its distributed pretty normal, thanks to Central Limit Theorem?