kubostat2018d p.2 :? bod size x and fertilization f change seed number? : a statistical model for this example? i response variable seed number : { i

Similar documents
k3 ( :07 ) 2 (A) k = 1 (B) k = 7 y x x 1 (k2)?? x y (A) GLM (k

講義のーと : データ解析のための統計モデリング. 第5回

kubostat2017c p (c) Poisson regression, a generalized linear model (GLM) : :

k2 ( :35 ) ( k2) (GLM) web web 1 :

kubostat2017b p.1 agenda I 2017 (b) probability distribution and maximum likelihood estimation :

講義のーと : データ解析のための統計モデリング. 第3回

一般化線形 (混合) モデル (2) - ロジスティック回帰と GLMM

kubostat2015e p.2 how to specify Poisson regression model, a GLM GLM how to specify model, a GLM GLM logistic probability distribution Poisson distrib

kubostat7f p GLM! logistic regression as usual? N? GLM GLM doesn t work! GLM!! probabilit distribution binomial distribution : : β + β x i link functi

(2/24) : 1. R R R

kubostat2017e p.1 I 2017 (e) GLM logistic regression : : :02 1 N y count data or

1 15 R Part : website:

1 Stata SEM LightStone 4 SEM 4.. Alan C. Acock, Discovering Structural Equation Modeling Using Stata, Revised Edition, Stata Press 3.

講義のーと : データ解析のための統計モデリング. 第2回

1 Stata SEM LightStone 3 2 SEM. 2., 2,. Alan C. Acock, Discovering Structural Equation Modeling Using Stata, Revised Edition, Stata Press.

(lm) lm AIC 2 / 1

Use R

Stata11 whitepapers mwp-037 regress - regress regress. regress mpg weight foreign Source SS df MS Number of obs = 74 F(

80 X 1, X 2,, X n ( λ ) λ P(X = x) = f (x; λ) = λx e λ, x = 0, 1, 2, x! l(λ) = n f (x i ; λ) = i=1 i=1 n λ x i e λ i=1 x i! = λ n i=1 x i e nλ n i=1 x

60 (W30)? 1. ( ) 2. ( ) web site URL ( :41 ) 1/ 77

y i OLS [0, 1] OLS x i = (1, x 1,i,, x k,i ) β = (β 0, β 1,, β k ) G ( x i β) 1 G i 1 π i π i P {y i = 1 x i } = G (

JOURNAL OF THE JAPANESE ASSOCIATION FOR PETROLEUM TECHNOLOGY VOL. 66, NO. 6 (Nov., 2001) (Received August 10, 2001; accepted November 9, 2001) Alterna

12/1 ( ) GLM, R MCMC, WinBUGS 12/2 ( ) WinBUGS WinBUGS 12/2 ( ) : 12/3 ( ) :? ( :51 ) 2/ 71

kubostat2017j p.2 CSV CSV (!) d2.csv d2.csv,, 286,0,A 85,0,B 378,1,A 148,1,B ( :27 ) 10/ 51 kubostat2017j (

最小2乗法

Microsoft Word - 計量研修テキスト_第5版).doc

第11回:線形回帰モデルのOLS推定

% 10%, 35%( 1029 ) p (a) 1 p 95% (b) 1 Std. Err. (c) p 40% 5% (d) p 1: STATA (1). prtesti One-sample test of pr

¥¤¥ó¥¿¡¼¥Í¥Ã¥È·×¬¤È¥Ç¡¼¥¿²òÀÏ Âè2²ó

/22 R MCMC R R MCMC? 3. Gibbs sampler : kubo/

9 8 7 (x-1.0)*(x-1.0) *(x-1.0) (a) f(a) (b) f(a) Figure 1: f(a) a =1.0 (1) a 1.0 f(1.0)

Microsoft Word - 計量研修テキスト_第5版).doc

こんにちは由美子です

untitled

yamadaiR(cEFA).pdf

R John Fox R R R Console library(rcmdr) Rcmdr R GUI Windows R R SDI *1 R Console R 1 2 Windows XP Windows * 2 R R Console R ˆ R

dvi

AR(1) y t = φy t 1 + ɛ t, ɛ t N(0, σ 2 ) 1. Mean of y t given y t 1, y t 2, E(y t y t 1, y t 2, ) = φy t 1 2. Variance of y t given y t 1, y t

<4D F736F F D20939D8C7689F090CD985F93C18EEA8D758B E646F63>

28

浜松医科大学紀要

第13回:交差項を含む回帰・弾力性の推定

2

25 II :30 16:00 (1),. Do not open this problem booklet until the start of the examination is announced. (2) 3.. Answer the following 3 proble

GLM PROC GLM y = Xβ + ε y X β ε ε σ 2 E[ε] = 0 var[ε] = σ 2 I σ 2 0 σ 2 =... 0 σ 2 σ 2 I ε σ 2 y E[y] =Xβ var[y] =σ 2 I PROC GLM

橡表紙参照.PDF

Microsoft Word - Win-Outlook.docx

/ 60 : 1. GLM? 2. A: (pwer functin) x y?

Microsoft Word - 計量研修テキスト_第5版).doc

こんにちは由美子です

24 Depth scaling of binocular stereopsis by observer s own movements

カテゴリ変数と独立性の検定

fx-9860G Manager PLUS_J

塗装深み感の要因解析

untitled

1 # include < stdio.h> 2 # include < string.h> 3 4 int main (){ 5 char str [222]; 6 scanf ("%s", str ); 7 int n= strlen ( str ); 8 for ( int i=n -2; i

R による共和分分析 1. 共和分分析を行う 1.1 パッケージ urca インスツールする 共和分分析をするために R のパッケージ urca をインスツールする パッケージとは通常の R には含まれていない 追加的な R のコマンドの集まりのようなものである R には追加的に 600 以上のパッ

/ 55 2 : : (GLM) 1. 1/23 ( )? GLM? (GLM ) 2.! 1/25 ( ) ffset (GLM )

,, Poisson 3 3. t t y,, y n Nµ, σ 2 y i µ + ɛ i ɛ i N0, σ 2 E[y i ] µ * i y i x i y i α + βx i + ɛ i ɛ i N0, σ 2, α, β *3 y i E[y i ] α + βx i


σ t σ t σt nikkei HP nikkei4csv H R nikkei4<-readcsv("h:=y=ynikkei4csv",header=t) (1) nikkei header=t nikkei4csv 4 4 nikkei nikkei4<-dataframe(n

untitled


Studies of Foot Form for Footwear Design (Part 9) : Characteristics of the Foot Form of Young and Elder Women Based on their Sizes of Ball Joint Girth

kubostat2018a p.1 統計モデリング入門 2018 (a) The main language of this class is 生物多様性学特論 Japanese Sorry An overview: Statistical Modeling 観測されたパターンを説明する統計モデル

統計モデリング入門 2018 (a) 生物多様性学特論 An overview: Statistical Modeling 観測されたパターンを説明する統計モデル 久保拓弥 (北海道大 環境科学) 統計モデリング入門 2018a 1

p.1/22

DAA09

ECCS. ECCS,. ( 2. Mac Do-file Editor. Mac Do-file Editor Windows Do-file Editor Top Do-file e

IPSJ SIG Technical Report Pitman-Yor 1 1 Pitman-Yor n-gram A proposal of the melody generation method using hierarchical pitman-yor language model Aki

Vol. 29, No. 2, (2008) FDR Introduction of FDR and Comparisons of Multiple Testing Procedures that Control It Shin-ichi Matsuda Department of

ブック

¥¤¥ó¥¿¡¼¥Í¥Ã¥È·×¬¤È¥Ç¡¼¥¿²òÀÏ Âè2²ó

こんにちは由美子です

4 OLS 4 OLS 4.1 nurseries dual c dual i = c + βnurseries i + ε i (1) 1. OLS Workfile Quick - Estimate Equation OK Equation specification dual c nurser

AN 100: ISPを使用するためのガイドライン

Duality in Bayesian prediction and its implication

紀要1444_大扉&目次_初.indd

山形大学紀要

1 環境統計学ぷらす 第 5 回 一般 ( 化 ) 線形混合モデル 高木俊 2013/11/21

3. ( 1 ) Linear Congruential Generator:LCG 6) (Mersenne Twister:MT ), L 1 ( 2 ) 4 4 G (i,j) < G > < G 2 > < G > 2 g (ij) i= L j= N

Rによる計量分析:データ解析と可視化 - 第3回 Rの基礎とデータ操作・管理

Page 1 of 6 B (The World of Mathematics) November 20, 2006 Final Exam 2006 Division: ID#: Name: 1. p, q, r (Let p, q, r are propositions. ) (10pts) (a

²¾ÁÛ¾õ¶·É¾²ÁË¡¤Î¤¿¤á¤Î¥Ñ¥Ã¥±¡¼¥¸DCchoice ¡Ê»ÃÄêÈÇ¡Ë

Motivation and Purpose There is no definition about whether seatbelt anchorage should be fixed or not. We tested the same test conditions except for t

第 55 回自動制御連合講演会 2012 年 11 月 17 日,18 日京都大学 1K403 ( ) Interpolation for the Gas Source Detection using the Parameter Estimation in a Sensor Network S. T

橡ボーダーライン.PDF

計量経済分析 2011 年度夏学期期末試験 担当 : 別所俊一郎 以下のすべてに答えなさい. 回答は日本語か英語でおこなうこと. 1. 次のそれぞれの記述が正しいかどうか判定し, 誤りである場合には理由, あるいはより適切な 記述はどのようなものかを述べなさい. (1) You have to wo

カルマンフィルターによるベータ推定( )

Read the following text messages. Study the names carefully. 次のメッセージを読みましょう 名前をしっかり覚えましょう Dear Jenny, Iʼm Kim Garcia. Iʼm your new classmate. These ar

10:30 12:00 P.G. vs vs vs 2

™…

1 kawaguchi p.1/81

201711grade2.pdf

LCC LCC INOUE, Gaku TANSEI, Kiyoteru KIDO, Motohiro IMAMURA, Takahiro LCC 7 LCC Ryanair 1 Ryanair Number of Passengers 2,000,000 1,800,000 1,

untitled

open / window / I / shall / the? something / want / drink / I / to the way / you / tell / the library / would / to / me


Stata 11 Stata ts (ARMA) ARCH/GARCH whitepaper mwp 3 mwp-083 arch ARCH 11 mwp-051 arch postestimation 27 mwp-056 arima ARMA 35 mwp-003 arima postestim

i5 Catalyst Case Instructions JP

Transcription:

kubostat2018d p.1 I 2018 (d) model selection and kubo@ees.hokudai.ac.jp http://goo.gl/76c4i 2018 06 25 : 2018 06 21 17:45 1 2 3 4 :? AIC : deviance model selection misunderstanding kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 1 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 2 / 44 number of parameters? 4 GLM 5 GLM : : 2012 05 18 http://goo.gl/ufq2 seed number 2 4 6 8 10 12 14 (A) k = 1 (B) k = 7 Too few parametes? 7 8 9 10 11 12 bod size x 2 4 6 8 10 12 14 Too man parameters? 7 8 9 10 11 12 bod size x How man parameters do ou need for the best prediction? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 3 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 4 / 44 :? :? k? 1. :? seed number 2 4 6 8 10 12 14 (A) k = 1 (B) k = 7 Too few parametes? 2 4 6 8 10 12 14 Too man parameters? 7 8 9 10 11 12 7 8 9 10 11 12 bod size x bod size x? number of parameters k? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 5 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 6 / 44

kubostat2018d p.2 :? bod size x and fertilization f change seed number? : a statistical model for this example? i response variable seed number : { i } explanator variable : bod size {x i } fertilization {f i } sample size f i C: T: control (f i = C): 50 sample (i {1, 2, 50}) fertilization (f i = T): 50 sample (i {51, 52, 100}) i x i probabilit distribution Poisson distribution : linear predictor : β 1 + β 2 x i + β 3 f i link function : log link function kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 7 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 8 / 44 :? 4 candidate models 4 : (A) constant λ :? 4 candidate models 4 : (B) f model λ i = exp(β 1 ) λ i = exp(β 1 + β 3 f i ) (log likelihood) (log likelihood) > loglik(glm( ~ 1, data = d, famil = poisson)) log Lik. 237.64 (df=1) > loglik(glm( ~ f, data = d, famil = poisson)) log Lik. 237.63 (df=2) kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 9 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 10 / 44 :? 4 candidate models 4 : (C) x model :? 4 candidate models 4 : (D) x + f model λ i = exp(β 1 + β 2 x i ) λ i = exp(β 1 + β 2 x i + β 3 f i ) (log likelihood) (log likelihood) > loglik(glm( ~ x, data = d, famil = poisson)) log Lik. 235.39 (df=2) > loglik(glm( ~ x + f, data = d, famil = poisson)) log Lik. 235.29 (df=3) kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 11 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 12 / 44

kubostat2018d p.3 :? k increases log L increases AIC : deviance (A) constant λ k = 1 237.6 (C) x model k = 2 235.4 (B) f model k = 2 237.6 fertilization Control (D) x + f model k = 3 235.3 Control fertilization 2. AIC : deviance badness of prediction : AIC kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 13 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 14 / 44 AIC output R glm() deviance : deviance AIC deviance D = 2 log L : deviance > glm( ~ x + f, data = d, famil = poisson) Call: glm(formula = ~ x + f, famil = poisson, data = d) Maximum log likelihood log L : goodness of fit Deviance D = 2 log L : Coefficients: (Intercept) x ft 1.2631 0.0801 0.0320 Degrees of Freedom: 99 Total (i.e. Null); Null Deviance: 89.5 Residual Deviance: 84.8 AIC: 477 97 Residual model k log L Deviance Residual 2 log L deviance constant λ 1 237.6 475.3 89.5 f 2 237.6 475.3 89.5 x 2 235.4 470.8 85.0 x + f 3 235.3 470.6 84.8 saturation 100 192.9 385.8 0.0 Residual Deviance? Null Deviance? AIC? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 15 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 16 / 44 AIC : deviance Null deviance, Residual deviance,... AIC : deviance badness of prediction : AIC = 2 log L + 2k Max deviance 475.3 470.8 constant λ x model Look for a model of the smallest AIC AIC Deviance 2 log L () 89.5 (Null Deviance) 85.0 (Residual Deviance) model k log L Deviance Residual 2 log L deviance AIC constant λ 1 237.6 475.3 89.5 477.3 f 2 237.6 475.3 89.5 479.3 x 2 235.4 470.8 85.0 474.8 x + f 3 235.3 470.6 84.8 476.6 saturation 100 192.9 385.8 0.0 585.8 Min deviance 385.8 saturation model AIC: A (or Akaike) information criterion kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 17 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 18 / 44

AIC : deviance (estimation)? constant λ ˆβ 1 = 2.04 β 1 = 2.08 parameter estimation kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 19 / 44 AIC : deviance Is it OK? Goodness of fit is evaluated b using the SAME data set...? constant λ ˆβ 1 = 2.04 log L () biased goodness of fit! kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 20 / 44 AIC : deviance : constant λ ˆβ 1 = 2.04 E(log L) β 1 = 2.08 ( ) (200 ) kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 21 / 44 AIC : deviance 1.95 2.00 2.05 2.10 2.15 2.20 140 135 130 125 120 115 110 log likelihood 1.95 2.00 2.05 2.10 2.15 2.20 140 135 130 125 120 115 110 1.95 2.00 2.05 2.10 2.15 2.20 140 135 130 125 120 115 110 β 1 β 1 β 1 (200 ) ( ) (A) (B) (A) (C) log L = 120.6 E(log L) = 122.9 ˆβ 1 = 2.04 β 1 = 2.08 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 22 / 44 AIC : deviance 1 2 2 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 23 / 44 3. likelihood ratio test kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 24 / 44 kubostat2018d p.4

kubostat2018d p.5 Although their procedures are similar... the are total different! AIC ( ) () AIC model selection totall different in their objectives kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 25 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 26 / 44 Objective? model selection : Look for a model of better prediction : rejection of null hpothesis kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 27 / 44 統計学的な検定 (NemanPearson framework) statistical test Null Alternative hpothesis hpothesis 帰無仮説対立仮説 glm( ~ 1) どうでもいい 興味ない VS glm( ~ x) 重要! これを主張したい! 非対称性 asmmetricit? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 28 / 44 統計学的な検定 (NemanPearson framework) statistical test Null Alternative hpothesis hpothesis 帰無仮説対立仮説 (if...) glm( ~ 1) test! reject 棄却 VS glm( ~ x) support 支持 非対称性 asmmetricit? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 29 / 44 統計学的な検定 (NemanPearson framework) statistical test (if...) Null hpothesis 帰無仮説 glm( ~ 1) test! NOT reject VS Alternative hpothesis 対立仮説 glm( ~ x) Sa Nothing!? 非対称性 asmmetricit? kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 30 / 44

kubostat2018d p.6 The same example, again test statistics D 1,2 i i seed number i D: deviance x i bod size x i neglect fertilization treatment (!) x model D 2 = 470.8 constant λ D 1 = 475.3 difference in deviance D 1,2 = D 1 D 2 = 4.51 4.5 likelihood ratio? log L 1 L = log L 1 log L 2 2 model k log L Deviance 2 log L constant λ 1 237.6 D 1 = 475.3 x 2 235.4 D 2 = 470.8 null hpothesis alternative hpothesis asmmetricit in test Null hpothesis is junk :... et we are focousing onl on null hpothesis kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 31 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 32 / 44 How to make null model objective : null hpothesis rejection Null hpothesis is included in Alt hpothesis this is a nested model ( ) observerd D 1,2 = 4.5 () ( ) () () { i } λ i alternative hpothesis : log λ i = β 1 + β 2 x i null hpothesis : log λ i = β 1 ( ) significant not significant is... (Reject ) (Not reject ) TRUE Tpe I error (no problem) NOT true (no problem) Tpe II error asmmetricit in test evaluating onl TpeI error : kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 33 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 34 / 44 generate D 1,2 distribution D 1,2 : Suppose null hpothesis is TRUE! bootstrap likelihood test How to generate D 1,2 under is TRUE?! constant λ x model D 1,2 ( ˆβ 1 = 2.06 ) D 1,2 D 1,2 D 1,2 > d$.rnd < rpois(100, lambda = mean(d$)) > fit1 < glm(.rnd ~ 1, data = d, famil = poisson) > fit2 < glm(.rnd ~ x, data = d, famil = poisson) > fit1$deviance fit2$deviance generation of random numbers virtual data rpois() ( ) fitting GLM to the virtual data glm() kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 35 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 36 / 44

kubostat2018d p.7 You must define rejection region in advance sa, 5%? 5%? NOT significant significant (5%) A random D 1,2 generator in R get.dd < function(d) # { n.sample < nrow(d) #.mean < mean(d$) # d$.rnd < rpois(n.sample, lambda =.mean) fit1 < glm(.rnd ~ 1, data = d, famil = poisson) fit2 < glm(.rnd ~ x, data = d, famil = poisson) fit1$deviance fit2$deviance # } pb < function(d, n.bootstrap) { replicate(n.bootstrap, get.dd(d)) } kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 37 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 38 / 44 Generated distribution of D 1,2 = D 1 D 2 Probabilit{ D 1,2 4.5} = 38 1000 = 0.038 observed D 1,2 D 1,2 = 4.5 constant λ x model D 1,2 (R code is in the next page) > source("pb.r") # reading "pb.r" text file > dd12 < pb(d, n.bootstrap = 1000) > hist(dd12, 100) # to plot histogram > abline(v = 4.5, lt = 2) > sum(dd12 >= 4.5) [1] 38 socalled P value is 0.038. kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 39 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 40 / 44 In this case, null hpothesis is rejected In case that P > 0.05...? So we can state that alternative hpothesis x model is better than constant λ. i i seed number i can be accepted. D: deviance x i bod size x i x model D 2 = 470.8 constant λ D 1 = 475.3 You can conclude NOTHING! You can NOT state that constant λ (Null hpothesis) is better λ asmmetricit in stattest : Null hpothesis is never accepted kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 41 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 42 / 44

kubostat2018d p.8 model selection misunderstanding model selection misunderstanding 4. model selection misunderstanding FAQ http://hosho.ees.hokudai.ac.jp/~kubo/ce/faqmodelselection.html kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 43 / 44 kubostat2018d (http://goo.gl/76c4i) 2018 (d) 2018 06 25 44 / 44