calculer la courbe gaussienne correspondant à une liste

Question

J'ai une liste de données comme ci-dessous. Je veux effectuer une régression non linéaire courbe de Gauss d'ajustement entre mediums et compte pour chaque élément de ma liste et le rapport moyenne et écart-type

mylist<- structure(list(A = structure(list(breaks = c(-10, -9, 
-8, -7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4), counts = c(1L, 
0L, 1L, 5L, 9L, 38L, 56L, 105L, 529L, 2858L, 17L, 2L, 0L, 2L), 
    density = c(0.000276014352746343, 0, 0.000276014352746343, 
    0.00138007176373171, 0.00248412917471709, 0.010488545404361, 
    0.0154568037537952, 0.028981507038366, 0.146011592602815, 
    0.788849020149048, 0.00469224399668783, 0.000552028705492686, 
    0, 0.000552028705492686), mids = c(-9.5, -8.5, -7.5, -6.5, 
    -5.5, -4.5, -3.5, -2.5, -1.5, -0.5, 0.5, 1.5, 2.5, 3.5), 
    xname = "x", equidist = TRUE), .Names = c("breaks", "counts", 
"density", "mids", "xname", "equidist"), class = "histogram"), 
    B = structure(list(breaks = c(-7, -6, -5, 
    -4, -3, -2, -1, 0), counts = c(2L, 0L, 6L, 2L, 2L, 1L, 3L 
    ), density = c(0.125, 0, 0.375, 0.125, 0.125, 0.0625, 0.1875 
    ), mids = c(-6.5, -5.5, -4.5, -3.5, -2.5, -1.5, -0.5), xname = "x", 
     equidist = TRUE), .Names = c("breaks", "counts", "density", 
    "mids", "xname", "equidist"), class = "histogram"), C = structure(list(
     breaks = c(-7, -6, -5, -4, -3, -2, -1, 0, 1), counts = c(2L, 
     2L, 4L, 5L, 14L, 22L, 110L, 3L), density = c(0., 
     0., 0.0246913580246914, 0.0308641975308642, 
     0.0864197530864197, 0.135802469135802, 0.679, 
     0.0185185185185185), mids = c(-6.5, -5.5, -4.5, -3.5, 
     -2.5, -1.5, -0.5, 0.5), xname = "x", equidist = TRUE), .Names = c("breaks", 
    "counts", "density", "mids", "xname", "equidist"), class = "histogram")), .Names = c("A", 
"B", "C")) 

J'ai lu ce Fitting a density curve to a histogram in R mais c'est comment adapter une courbe à un histogramme. ce que je veux est des valeurs de meilleur ajustement »

« Mean » « SD »

Si je PRISM de le faire, je devrais obtenir les résultats suivants A

Mids Counts 
-9.5 1 
-8.5 0 
-7.5 1 
-6.5 5 
-5.5 9 
-4.5 38 
-3.5 56 
-2.5 105 
-1.5 529 
-0.5 2858 
0.5  17 
1.5  2 
2.5  0 
3.5  2 

exécution non linéaire ajustement de courbe gaussienne régression, je me

"Best-fit values" 
"  Amplitude" 3537 
"  Mean"  -0.751 
"  SD"   0.3842 

pour le second ensemble B

Mids Counts 
-6.5 2 
-5.5 0 
-4.5 6 
-3.5 2 
-2.5 2 
-1.5 1 
-0.5 3 



"Best-fit values" 
"  Amplitude" 7.672 
"  Mean"   -4.2 
"  SD"   0.4275 

et pour la troisième

Mids Counts 
-6.5 2 
-5.5 2 
-4.5 4 
-3.5 5 
-2.5 14 
-1.5 22 
-0.5 110 
0.5  3 

Je reçois ce

"Best-fit values" 
"  Amplitude" 120.7 
"  Mean"  -0.6893 
"  SD"  0.4397 

Answer 1

Afin de convertir l'histogramme retour à l'estimation de la moyenne et l'écart type. Convertissez d'abord les résultats du nombre de casiers par bin. Ce sera une approximation des données d'origine.

Sur la base de votre exemple ci-dessus:

#extract the mid points and create list of simulated data 
simdata<-lapply(mylist, function(x){rep(x$mids, x$counts)}) 
#if the original data were integers then this may give a better estimate 
#simdata<-lapply(mylist, function(x){rep(x$breaks[-1], x$counts)}) 

#find the mean and sd of simulated data 
means<-lapply(simdata, mean) 
sds<-lapply(simdata, sd) 
#or use sapply in the above 2 lines depending on future process needs 

Si vos données entiers a ensuite été en utilisant les pauses que les bacs fourniront une meilleure estimation. Selon la fonction de l'histogramme (c.-à-d. Right = TRUE/FALSE), les résultats peuvent être décalés d'une unité.

Modifier

Je pensais que cela allait être facile. J'ai examiné la vidéo, les données d'échantillon présentées ont été:

mids<-seq(-7, 7) 
counts<-c(7, 1, 2, 2, 2, 5, 217, 70, 18, 0, 2, 1, 2, 0, 1) 
simdata<-rep(mids, counts) 

Les résultats vidéo étaient moyenne = -0,7359 et sd = 0,4571. La solution que j'ai trouvé fourni les résultats les plus proches utilisait le « fitdistrplus » package:

fitdist(simdata, "norm", "mge") 

Utilisation de la « maximisation de l'estimation de qualité d'ajustement » entraîné moyenne = -0,7597280 et sd = 0,8320465.
À ce stade, la méthode ci-dessus fournit une estimation proche, mais ne correspond pas exactement. Je ne sais pas quelle technique a été utilisée pour calculer l'ajustement de la vidéo.

Edit # 2

Les solutions ci-dessus impliqués recréer les données d'origine et de montage qui en utilisant soit la moyenne/SD ou en utilisant le package fitdistrplus. Cette tentative est une tentative d'effectuer un ajustement des moindres carrés en utilisant la distribution gaussienne.

simdata<-lapply(mylist, function(x){rep(x$mids, x$counts)}) 
means<-sapply(simdata, mean) 
sds<-sapply(simdata, sd) 

#Data from video 
#mids<-seq(-7, 7) 
#counts<-c(7, 1, 2, 2, 2, 5, 217, 70, 18, 0, 2, 1, 2, 0, 1) 

#make list of the bins and distribution in each bin 
mids<-lapply(mylist, function(x){x$mids}) 
dis<-lapply(mylist, function(x) {x$counts/sum(x$counts)}) 

#function to perform the least square fit 
nnorm<-function(values, mids, dis) { 
    means<-values[1] 
    sds<-values[2] 
    #print(paste(means, sds)) 
    #calculate out the Gaussian distribution for each bin 
    modeld<-dnorm(mids, means, sds) 
    #sum of the squares 
    diff<-sum((modeld-dis)^2) 
    diff 
} 

#use optim function with the mean and sd as initial guesses 
#find the mininium with the mean and SD as fit parameters 
lapply(1:3, function(i) {optim(c(means[[i]], sds[[i]]), nnorm, mids=mids[[i]], dis=dis[[i]])}) 

Cette solution fournit une réponse plus proche aux résultats de PRISM, mais pas toujours la même chose. Voici une comparaison de toutes les 4 solutions. À partir de la table, l'ajustement des moindres carrés (celui juste au-dessus) fournit l'approximation la plus proche. Peut-être que peaufiner la fonction dnorm des points milieu pourrait aider. Mais les données du cas B sont les plus éloignées d'une distribution normale, mais le logiciel PRISM génère toujours un petit écart-type, alors que les autres méthodes sont similaires. Il est possible que le logiciel PRISM effectue un certain type de filtrage de données pour supprimer les valeurs aberrantes avant l'ajustement.