0

I have a list data like below. I want to perform nonlinear regression Gaussian curve fitting between mids and counts for each element of my list and report mean and standard deviation 

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.0123456790123457, 
        0.0123456790123457, 0.0246913580246914, 0.0308641975308642, 
        0.0864197530864197, 0.135802469135802, 0.679012345679012, 
        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"))

I have read this Fitting a density curve to a histogram in R but this is how to fit a curve to a histogram. what I want is Best-fit values"

" Mean" " SD"

If I use PRISM to do it, I should get the following results for 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

performing nonlinear regression Gaussian curve fitting , I get 

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

for the second set 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

and for the third one 

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

I get this 

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

1 に答える 1

1

ヒストグラムを平均値と標準偏差の推定値に変換するため。まず、ビン カウントの結果をビンに変換します。これは、元のデータの近似値になります。

上記の例に基づいて：

#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

データが整数の場合、ブレークをビンとして使用すると、より適切な推定値が得られます。ヒストグラムの関数 (つまり、右 = TRUE/FALSE) によっては、結果が 1 つシフトする場合があります。

編集

これなら簡単にいけると思いました。ビデオを確認したところ、表示されたサンプル データは次のとおりです。

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)

ビデオの結果は、平均 = -0.7359 および sd = 0.4571 でした。最も近い結果が得られた解決策は、「fitdistrplus」パッケージを使用することでした。

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

「適合度推定の最大化」を使用すると、平均 = -0.7597280 および sd = 0.8320465 が得られました。
この時点で、上記の方法は正確な見積もりを提供しますが、正確には一致しません。ビデオからフィットを計算するためにどのような手法が使用されたかはわかりません。

編集＃2

上記の解決策には、元のデータを再作成し、平均/標準偏差または fitdistrplus パッケージを使用してフィッティングすることが含まれていました。この試行は、ガウス分布を使用して最小二乗近似を実行する試みです。

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]])})

このソリューションは、PRISM の結果により近い回答を提供しますが、それでも同じではありません。以下は、4 つのソリューションすべての比較です。 ここに画像の説明を入力

表から、最小二乗適合 (すぐ上のもの) が最も近い近似値を提供します。中間点の dnorm 関数を微調整すると役立つ場合があります。ただし、ケース B のデータは正規分布から最も離れていますが、PRISM ソフトウェアは依然として小さな標準偏差を生成しますが、他の方法は同様です。PRISM ソフトウェアは、当てはめの前に外れ値を除去するために何らかのデータ フィルタリングを実行する可能性があります。

于 2016-06-28T15:42:56.147 に答える