Excelによるパラメータ推定

このQ&Aのポイント
  • Excelを使用して、ハフモデルの魅力度を求め、パラメータの推定を行おうとしていますが、魅力度とパラメータが未知のため、推定方法について困っています。
  • ハフモデルの魅力度は、売場面積と人口の重みによって求められますが、これらの値は既知です。
  • しかし、パラメータが未知のため、推定方法がわからず混乱しています。どのように推定すればよいでしょうか?
回答を見る
  • ベストアンサー

Excelによるパラメータ推定

http://aoki2.si.gunma-u.ac.jp/Hanasi/StatTalk/solver.html このページを参考に、ハフモデルの魅力度を求め、パラメータの推定を行おうとしているのですが、 魅力度は未知であり、またパラメータも未知なのでどのようにして推定をしていいかわかりません。 魅力度の式は  Wj=Sj^θ1×Σ(Oi/Tij)^θ2 Wj;魅力度 Sj;売場面積 Oi;所要時間に関する人口の重み Tij;iからjへの所要時間 で、 既知な値はSj、Oi、Tij  です。 頭が混乱して質問の意図が明確になってないかもしれません・・・ よろしくお願いします。

  • jou-s
  • お礼率72% (72/100)

質問者が選んだベストアンサー

  • ベストアンサー
  • at9_am
  • ベストアンサー率40% (1540/3760)
回答No.2

#1です。 手順としては、以下のようになります。 a) アドインのソルバーをインストールしておく。 b) ln(Wj)、ln(Sj)、ln(Oi/Tij) を入れた票をエクセルで作る。 c) θの初期値として適当な数値(0など)を入れる。 d) 残差 e = ln(Wj) - {θ1 ln(Sj) + θ2 ln(Oi/Tij)} の二乗和を計算する。二乗和は sumsq 関数を使うと簡単に求められる。 e) 残差二乗和を最小にするようなパラメータを計算させる。ここについてはHPの手順7~9と同じなので割愛。 標準偏差は、求められたパラメータを用いた残差から計算するのですが、線形化された関数の最小自乗推定を行うのであれば、同じくアドインの分析ツールにある回帰分析の方が使い勝手が良いでしょう。こちらはyに被説明変数、xに説明変数を指定すれば標準偏差、t統計量、p値などが手軽に得られます。 ところで。 前回書き忘れましたが、この式は定数項が必要なのではありませんか? つまり魅力度の式は  Wj= θ0 × Sj^θ1×Σ(Oi/Tij)^θ2 ではありませんか?

jou-s
質問者

お礼

訂正がありました。 もともとのデータは売場面積、従業員数、店舗数、Oi*Tij の4つだけで、 魅力度の値とパラメータの値は分析により出力されたものだと思います。 (過去のデータなのでどうやって分析したのかわからない)

jou-s
質問者

補足

お返事ありがとうございます。 今回推定している魅力度は過去にある論文のものなのですが、 魅力度推定のエクセルデータで 店舗数、従業員、売場面積、Oi*Tij、魅力度、パラメータ(e1~e3)の数値データが打ち込まれており、 魅力度の式が  =売場面積^e1*(従業員数/店舗数)^e2*(Oi*Tij)^e3 に指定されていました。 パラメータe2の値が0になっていたので論文では省略されていたのかもしれません。 このデータからソルバー機能でパラメータを求める際に、 標準偏差などの値がなかったため困惑していました。 まわりくどくて申し訳ありません。。。

その他の回答 (1)

  • at9_am
  • ベストアンサー率40% (1540/3760)
回答No.1

与式に対数をとると log(Wj) = θ1 log(Sj) + Σθ2i (Oi/Tij) となります。 数値の対数変換はエクセルでもすぐにできます。 初期時点でのパラメータは、例えば全て0でも構いません。最小自乗法の場合、極小値が一つしかないことが保証されていますから、時間の差こそあれ、正しい答に到達します。

jou-s
質問者

補足

回答ありがとうございます。 式を対数変換してからパラメータ推定をするのですね。 標準偏差はどのようにしたらよいのでしょうか。 よろしければエクセルでの具体的な操作を教えていただけたら幸いです。

関連するQ&A

  • 統計学、最尤推定について

    統計学、最尤推定について こんにちは。早速ですが、質問させてください。 ある分布、ここでは例としてワイブル分布W(m,η)のパラメータについて考えます。 ここに、mはワイブルパラメータ、ηは尺度パラメータです。 ワイブル分布に従うデータx1,...,xnからパラメータを最尤推定します。 このとき、例えば、m=2が真の値だとして、その値を既知としたもとで最尤推定 された\hat{η}(m)と、その値を未知としたもとで最尤推定された\hat{m}の共分散を 求めたいがどうすれば良いかと言うのが質問です。 分かりづらいシチュエーションですが、ある証明の中で個の共分散を使いたいのです。 よろしくお願いいたします。

  • 正規分布への適合度検定について

    統計ソフト(R)において、正規分布への適合度検定をするための、カイ二乗検定を用いた関数(normaldist←青木先生という方のサイト上にありました。http://aoki2.si.gunma-u.ac.jp/R/normaldist.html)を使う際に、度数分布表における階級の分け方(どれくらい細かく階級をとるか)を変えると、検定結果のp値が変化するのですが、階級の決まった分け方というものはあるのでしょうか?よろしくお願いします。

  • 仮説検定のキモを教えてください。

    例えば下記のサイトの例で http://aoki2.si.gunma-u.ac.jp/lecture/Cross/cross.html 「例題では,自由度 6 の χ2 分布において,Pr{χ2≧12.59}=0.05 であるから,  P=Pr{χ2≧13.713}<0.05 である(正確な有意確率:P=0.03301)」 とありますが、2変数の独立性を無くすために 表1の観測値の一部もっと大きくすると、 ⇒期待値との差が広がる ⇒χ2 0の値は12.59よりも大きくなる。 ⇒Pr{χ2≧12.59より大きい値}<0.05 のまま ⇒帰無仮説は棄却される。 で、理解不能になりました。 どなたか教えてください。

  • 寿命予測について

    物品の経過年数別の撤去率もしくは生存率から寿命を予測する場合、 1.経過年数別の撤去率に成長曲線(ゴンペルツ曲線)をあてはめて寿命予測。 2.経過年数別の生存率にワイブル分布もしくは指数分布をあてはめて寿命予測。 以上2つの方法を比較した場合、どのような長所、短所があるのでしょうか? ヒトの生命表の老年部分の予測には、ゴンペルツメイカム曲線が使われ、 工学的な故障率の推定にはワイブル分布がよく使われているようですが。 http://aoki2.si.gunma-u.ac.jp/lecture/ 上記URLや「寿命の数理」古川俊之 著等を参照しているのですが 良く分かりません。

  • OPENCVを使ったエピポーラ幾何からのカメラ運動推定

    OPENCVでUSBカメラを動かしながら撮影した画像列から、 カメラがどのように動いたかを推定したいと考えています。 物体の3次元復元ではなく、カメラの6自由度の運動パラメータの復元がしたいのです。 FindFundamentalMat()を使って、2枚の画像間の対応点から 基礎行列(F行列)を求めることは出来ました。 しかし、このF行列から運動パラメータを復元する方法が分かりません。 参考書やウェブで調べても、F行列の求め方や、 F行列が回転・並進の運動パラメータと未知の内部パラメータの情報を含む、 といった説明はあるのですが、 そこから回転と並進の運動パラメータを復元する説明がなく、 「F行列からRとtが復元できる」と簡単に書かれているばかりです。 具体的にどのような操作をすれば、F行列から回転Rと並進tが求まるのでしょうか? また、最終的に得たい情報は「各軸回りの回転角は何度か?」と 「各軸方向にどれだけ動いたか?」であるわけですが、 F行列から求まるRはどのような形になるのでしょうか? 参考書ではRは『3×3の行列』だとしているのですが、 R= |1 0  0 | |cosY 0 sinY||cosZ -sinZ 0| |0 cosX -sinX| |0  1  0 ||sinZ cosZ 0| |0 sinX cosX| |-sinY 0 cosY||0  0  1| だともされています。 単純に、この3つの行列の掛け算結果が3×3行列の形で求まるのでしょうか? できれば、実例をあげて説明していただくか、 そのようなものが掲載されている書籍かサイト(できれば日本語で)を 教えていただけますでしょうか。 どうかよろしくお願いします。 なお、参考書として以下を利用しています。 「コンピュータビジョン-視覚の幾何学-」佐藤淳著 コロナ社刊 1999年

  • カイ二乗の自由度

    実験による測定値に対しカイ二乗を用いる際の自由度についての質問です。 xを指定した上で、yを測定するという作業を異なるxについてn回行い、n組の(x_n, y_n)を得ます。 このyは y = f(x, A, B, C)とモデル化された関数f()で書くことができ(x, A, B, Cは全て独立)、A, B, Cは未知のパラメーターです。 最終的にこのA, B, Cを推定することを目標としますが、ここで、nがあまり多く得られないので、 最小二乗などのあてはめにより最適なA, B, Cを求めることは諦め、適当な棄却値(%)を決めた上で、 カイ二乗によりA, B, Cが取り得る値の範囲を探すことを考えます。 実際にやることは、考えられる有限(i, j, k個)のA_i, B_j, C_k全ての組み合わせに対しカイ二乗を計算し、 その自由度での棄却値に対応するカイ二乗値より小さいカイ二乗値を示す点(A_i, B_j, C_k)を探します。 ここで質問ですが、求めるパラメータの数は3なのですが、この場合の自由度もn-3として良いのでしょうか? それとも、各点の計算ではA, B, Cを既知のものとして与えてしまうわけですので、 パラメータ数3を引かずに自由度nとするべきでしょうか? よろしくお願い致します。

  • 関数のパラメータ推定について

    http://bekkoame.okwave.jp/qa2595561.htmlで質問させて頂いた者です。 y = γ / (1 + exp(α + β * x))という関数において,パラメータγ,α,βを推定したいのですが,初期値の設定がうまくいかなくて困っております。要するに初期値の探し方を教えて欲しいのです。 この関数からγ / y -1 = exp(α - β * x)という式が導けます。ここで両辺の自然対数をとると, log(γ / y - 1) = α - β * x となります。ここでY = log(γ / y - 1)とするとY = α - β * xとなって,これは線形モデルなので最小二乗法によって簡単にαとβの値が求められます。今回の場合だと得られているデータは, x=[54,57,60,63,66,69,72] y=[0,1,2,5,8,9,10] なので,とりあえずγ=11と仮定して計算すると結局は Y = [*, 2.30, 1.50, 0.18, -0.98, -1.50, -2.30] となります(*はうまく計算できないので無視)。最終的にYのデータからY(ハット) = -0.949 * x + 3.1883という式が得られます。以前の質問でも補足しましたが,Rというデータ解析用のソフトウェアを使うと, nls(y ~ gamma / ( 1 + exp(alpha + beta * x)),start = c(gamma = 11,alpha = 3.1883,beta = -0.949)) とするのですが,これは要するにy ~ gamma / ( 1 + exp(alpha + beta * x))の部分がモデル式,start = c(gamma = 11,alpha = 3.1883,beta = -0.949)の部分が初期値を設定する部分です。yとxという変数には上に記したように得られたデータがベクトルとして格納されています。 以上,説明が長くなってしまいましたが,今回のようにx=[54,57,60,63,66,69,72]となっていると,初期値が悪いとのエラーが出て計算ができません。これをx=[1,2,3,4,5,6,7]としてやるとうまくいきます。だから平行移動させるために前回のような質問をしたのですが、、、(ちょっと文字数に制限があるので意味不明だったら補足します)

  • 回帰係数の有意性の検定

    線形回帰式 y=b0+Σ[i=1~p]bi・xi y:被説明変数 xi:説明変数 bi:係数 があり、N組のデータからbiの最尤推定量 bi~ を求めるとします。 この bi~ の有意性を検定する場合、bi~ が最尤推定量の漸近正規性から、以下の正規分布に従っているとみなします。 N(bi★,σ^2) (bi★は真のパラメータ、分散σ^2は未知) 次に以下の仮説検定を行います。 帰無仮説 H0:bi★=0 対立仮説 H1:bi★≠0 ここでは bi~ の分散が未知なので、bi~ の標本分散 s^2 を用いて 以下のt統計量により検定を行います。 t=bi~/(s/√N) ~ 自由度N-1のt分布 以上のような解説を聞いたのですが 最尤推定量 bi~ はN個のデータから1回推定を行っただけなので 標本は1個しかありませんよね? だったら、標本分散s^2って計算できない気がするのですが… 私の理解のどこがおかしいでしょうか?

  • 対数尤度からのパラメター推定

    式(1)の対数尤度から未知数を推定するためにパラメターxで偏微分すると式(2)になるようなのですが、式(1)から(2)への過程を解説していただけますでしょうか。 logL(x|y,s^2) = Sigma(i=1 to n)logP([y-fx]/s) (1) d log(x|y,s^2) ----------- = Sigma(i=1 to n) (Aij/s)W(v/s) (2) d x Aij はヤコビアン係数W df/dx(偏微分)、vは残差  関数W(z)は W(z) = -dlogP(z)/dz よろしく御願いいたします。

  • MCMCを用いたパラメータ推定[Matlab]

    今、MCMCを用いてパラメータ推定のためのプログラムを考えています。Matlabを用いてコードを書いている最中なのですが、IfとEndが合わないという忠告が何個も出てしまいます。下記にコードを載せましたので、何が間違っているか教えていただけませんでしょうか。 求めたいパラメータは、ポアッソン分布に従うと仮定し,尤度を使ってデータの当てはまりを見ています。 Matlab初心者なので間違いも多いと思いますが、どうぞよろしくお願いいたします。 logpt =0; %空値 theta =0; %空値 theta0 = ceil(100.*rand()); % i = initial value of MCMC 100 lambda = 5; % true mean value of dispersal updatewidth = 0.1; %sigma update width of MCMC nsample = 100; %number of iterations for MCMC nsample 100000 nthin = 1; %thinning number of MCMC kannbiki, sukashi individual = 100; %individual numbers 1000 xsize = 10; %space size of matrix ysize = 10; dat = zeros(xsize, ysize); %space xinit = 5; %initial position of population yinit = 5; %making dammy data for j = 1:individual lam(1) = ceil(100.*rand()); % decide a initial number P = poissrnd(lambda); tx = xinit; ty = yinit; i = 1; while i <= P d = ceil(4.*rand()); if d==1 && ty -1 >0 ty = ty-1; i = i+1; elseif 1==2 && tx +1 <= xsize tx = tx+1; i = i+1; elseif d==3 && ty+1 <= ysize ty = ty+1; i = i+1; elseif d==4 && tx-1 > 0 tx = tx-1; i = i+1; end; end; dat(tx,ty) = dat(tx,ty) + (1); end; %MCMC start for n = 1:nsample for o = 1:nthin %(ii) simulation with the previous step dispersal value sp0 = zeros(xsize, ysize); ind0 = zeros(2,individual); for m = 1:individual r = poissrnd(theta0); tx = xinit; ty = yinit; i = 1; while i <= r d = ceil(4.*rand()); if d==1 && ty-1 > 0 ty = ty-1; i = i+1; elseif d==2 && tx+1 <= xsize tx = tx+1; i = i+1; elseif d==3 && ty+1 <= ysize ty = ty+1; i = i+1; elseif d==4 && tx-1 > 0 tx = tx-1; i = i+1; end; end; sp0(tx,ty) = sp0(tx,ty) + 1; ind0(1,m) = tx; ind0(2,m) = ty; end; %calculating likelihood LOGP = 0; for q = 1:xsize for s = 1:ysize if sp0(q,s) == 0 LOGP = LOGP + log(poisspdf(0.001,dat(q,s))); else LOGP = LOGP + log(poisspdf(sp0(q,s),dat(q,s))); end; end; end; %(iii) making a proposal value tmp = updatewidth * randn(size(theta0)); %b= updatewidth = 0.1; sigma update width of MCMC while tmp < 0 tmp = updatewidth * randn(size(theta0)); end; %(iv) making proposal distances of each individuals theta1 = poissrnd(tmp); %(v) simulation with the proposal value for p = 1:individual tx = xinit; ty = yinit; i = 1; while i <= theta1 d = ceil(4.*rand()); if d==1 && ty-1>0 ty = ty-1; i = i+1; elseif d==2 && tx+1 <= xsize tx = tx+1; i = i+1; elseif d==3 && ty+1 <= ysize ty = ty+1; i = i+1; elseif d==4 && tx-1 > 0 tx = tx-1; i = i+1; end; end; %calculating likelihood prev = sp0(ind0(1,p),ind0(2,p)); dat_prev = dat(ind0(1,p),ind0(2,p)); prop = sp0(tx,ty); dat_prop = dat(tx,ty); if prop == 0 LOGP0 = log(poisspdf(prev,dat_prev)) + log(poisspdf(0.001,dat_prop)); else LOGP0 = log(poisspdf(prev,dat_prev)) + log(poisspdf(prop,dat_prop)); end; if prev == 1 LOGP1 = log(poisspdf(0.001,dat_prev)) + log(poisspdf(prop+1,dat_prop)); else LOGP1 = log(poisspdf(prev-1,dat_prev)) + log(poisspdf((prop+1),dat_prop)); end; if rand < exp(LOGP1 - LOGP0) r = theta1; LOGP0 = LOGP1; end; end; %(vi) calculating the mean of theta theta0 = mean(r); end; theta = theta0; LOGPt = LOGP; end;