2016-06-28 6 views
0

私は以下のようなリストデータを持っています。私はこの Fitting a density curve to a histogram in R を読んだことがあるが、これはフィットする方法ですリストにガウス曲線近似を計算する

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")) 

意味私のリストおよびレポートの各要素に対して 中音域カウント間の非線形回帰近似ガウス曲線を実行すると標準偏差ヒストグラムにカーブします。私が欲しいのは、私はそれを行うにはPRISMを使用している場合は、私は非線形を実行

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 

ための次のような結果 を取得する必要があります 『 『SD』

平均「

』ベストフィット値であります回帰ガウス曲線フィッティング、私は

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

は、第二セット 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 

と第三1

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 

のために私は戻って平均値と標準偏差の推定値にヒストグラムを変換するために、この

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

推定平均と標準偏差/分散を求めているなら、これは最尤法で達成できると思います。ベースRと 'maxLik'パッケージに' mle'関数があります。この例では、ミッドとカウントではなく、生データを使用する必要があります。 'mle'の最初の例は、あなたが望むものと類似したものでなければなりません。 – lmo

+0

私は現時点でビデオを見ることができませんが、私ができる時は数時間後にそれを見せます。ビニングされたデータから推定すると有用な情報が失われているようです。これは、あなたがそのような小さな標本サイズを持っていることを考えると、特にそうです。 – lmo

+0

@lmo大丈夫ですが、実際にはサンプルのサイズは1000ほど遥かに高くなっています。その場合は問題にならないと思います – nik

答えて

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 

。ヒストグラムの関数(right = 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でした。 =平均の結果-0.7597280とSD = 0.8320465

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

「最大適合度推定」を使用する:私が見つけた解決策は、最も近い結果が「fitdistrplus」パッケージを使用していた提供しました。
この時点で、上記の方法は近い見積もりを提供しますが、正確には一致しません。ビデオからのフィットを計算するためにどのようなテクニックが使われたのか分かりません。

編集#2

上記の解決策は、元のデータとフィッティング平均/ SD又は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つのソリューションのすべてを比較しています。 enter image description here

表から、最小二乗フィット(直上のもの)が最も近似しています。多分中点を微調整すると、dnorm関数が役に立ちます。しかし、ケースBのデータは正規分布から最も離れているが、PRISMソフトウェアは小さな標準偏差を生成しているが、他の方法は類似している。 PRISMソフトウェアは、適合の前に異常値を除去するために何らかのタイプのデータフィルタリングを実行することが可能である。

+0

これを行うことで、1つは非線形回帰ガウス曲線適合を適用するのですか? – nik

+0

こんにちは、上記のデータを確認し、PRISMを使用して非線形回帰ガウス曲線近似を作成し、平均と標準偏差を求めました。それが同じであるかどうか見ていただけますか? – nik

+0

値が一致しません。 PRISMソフトウェアがどのように機能しているかわかりません。フィット感のために尾を切り取ったり平滑にしたりするかもしれません。あなたのケースBはあまり正常ではありませんがPRISMは<0.5の標準偏差を生成しています – Dave2e

関連する問題