Elección del centro de las franjas del histograma para el ajuste

Question

Tengo una distribución bimodal, y si se traza con Mathematica se ve así:

Ahora, el valor más bajo de los datos reales es 8196 y 690720, pero como se ve en el gráfico, Mathematica deja que el rango de datos vaya de 0 a 744572. ¿Está Mathematica eligiendo un mal rango de datos para el histograma?

¿Cuál es, en general, una buena opción para definir el centro de los bins y el rango de datos, de manera que pueda ajustar una distribución a través del histograma?

Mi enfoque sería:

(1) Calculate (bin width) = (Max-Min)/(number of bins)
    [I'm aware that there are different rules how to choose the optimal number of
    bins depending on the underlying distribution, let's just assume this is 12]

(2) then I have 12 equal bins, starting from 8196 and ending at 690720, each having a width of 56877

(3) The first bin goes then from 8196 to 65073=8196+56877 and so on

(4) As the center of the bin I define the middle between 8196 and 65073 which is 36624.5 and I position my first bin there.

(5) Then I get 12 data pairs of bin center position and number of observations and I can fit a bimodal distribution through it

¿Estoy cometiendo un error si hago eso, o cuál es el razonamiento detrás de la elección de Mathematica de que el rango del histograma exceda el rango real de los datos?

Edición: He subido los datos en bruto aquí: datos en bruto

Edit2: Para aclarar la misteriosa frecuencia de 17 que fue señalada por Stephan: Los datos son un mapa de fotoluminiscencia confocal donde un láser escanea un emisor y se ve así: El láser escanea fila por fila, por lo que en medio de cada fila se enciende el emisor, lo que explica la frecuencia de 17 cuando se trazan los datos brutos, ya que se originan en una sola lista.

Answer 1

Como @Nick Cox dice , ajuste su distribución directamente a los datos. No ponga primero los datos en un histograma. ¿Por qué querría hacerlo?

En su lugar, ajuste una densidad de núcleo estándar. Usaré R, porque lo conozco mejor, pero asumo que Mathematica tiene funcionalidades similares. (Si no es así, te recomiendo que aprendas R.) A continuación hay un código que ajustará dicha densidad a tus datos y extraerá el $x$ para el segundo pico.

Para mayor claridad, podemos evaluar el grado de seguridad de esta coordenada haciendo un bootstrap. También estoy trazando un cuantil del 95% de la prueba de arranque. Observe cómo esto es ligeramente asimétrico.

dataset <- unlist(read.table("https://files.fm/down.php?i=qucxqxgw"))

foo <- density(dataset)
max.index <- which(foo$x>4e5)[which.max(foo$y[foo$x>4e5])]

plot(foo)
points(foo$x[max.index],foo$y[max.index],pch=19,col="red")
text(foo$x[max.index],foo$y[max.index],round(foo$x[max.index]),pos=3,col="red")

library(boot)
bootstrap <- boot(dataset,statistic=function(dataset,index){
    foo <- density(dataset[index])
    max.index <- which(foo$x>4e5)[which.max(foo$y[foo$x>4e5])]
    foo$x[max.index]
}, R=1e3)

lines(quantile(bootstrap$t,c(0.025,0.975)),rep(foo$y[max.index],2),col="red",lwd=2)

Si quiere el anchura del segundo pico, se puede extraer de la densidad (y hacer un bootstrap) después de haber decidido cómo se define un pico (algo más del 95% del valor del pico, o un desplazamiento fijo, o algo más).

(Sí, en principio podríamos corregir el hecho de que sus datos parecen ser todos no negativos, mientras que la estimación de la densidad es negativa. En la práctica, dado que sólo te interesa el segundo pico, no veo realmente el sentido).

Sin embargo...

Aquí hay un gráfico de sus datos originales:

plot(dataset,type="o")

Esto parece extrañamente regular. Así que, después de jugar un poco con el frequency encontramos el siguiente diagrama de estaciones:

library(forecast)
seasonplot(ts(dataset,frequency=17))

A menos que hayas hecho una ordenación muy extraña de tus datos brutos, tus datos son realmente estacionales con un periodo de 17. Por lo tanto, me pregunto si encontrar la ubicación de su segundo modo en tales datos es realmente lo que usted quiere hacer en absoluto.

Answer 2

Puedes usar Mathematica para dibujar un SmoothHistogram. Aquí uso los datos que me proporcionaste en una lista aplanada.

data = Flatten[Import["raw.txt", "TSV"]]

Histograma[datos]

SmoothHistogram[data]

Ahora suponga que quiere ajustar a una mezcla de dos distribuciones, y obtener los parámetros y también la mezcla, primero defina una mezcla de distribución:

distMix = MixtureDistribution[{p, 1 - p}, {NormalDistribution[a, b], NormalDistribution[c, d]}]

Y luego, obtener los parámetros:

params = FindDistributionParameters[data, distMix]

{p -> 0.382944, a -> 516841., b -> 101764., c -> 124503., d -> 81260.4}