Bayesian Variable Selection Method for Modeling Dose-Response - - PowerPoint PPT Presentation

Bayesian Variable Selection Method for Modeling Dose-Response Microarray Data Under Simple Order Restrictions

Bayes2013, Rotterdam

Martin Otava

I-Biostat, Universiteit Hasselt

22.05.2013

Martin Otava (I-Biostat, UHasselt) BVS 1 / 20

Research team

Hasselt University (Belgium):

Martin Otava Ziv Shkedy

Durham University (UK):

Adetayo Kasim

Imperial College London (UK):

Bernet Kato

Zoetis (Belgium):

Dan Lin

Janssen Pharmaceutica (Belgium):

Luc Bijnens Hinrich W.H. G¨

• hlmann

Willem Talloen

Martin Otava (I-Biostat, UHasselt) BVS 2 / 20

Model formulation

Dose-response modeling

Increasing dose of therapeutical compound. Variety of possible responses:

Toxicity. Inhibition or stimulation. Gene expression level.

Goal:

Determine if there is any relationship. If so, what is the shape of the profile. Select threshold doses (e.g. MED).

• Dose

Gene expression

1 2 3 7.0 7.2 7.4 7.6 7.8

• Dose

Gene expression

1 2 3 6.0 6.1 6.2 6.3 6.4 6.5

Martin Otava (I-Biostat, UHasselt) BVS 3 / 20

Model formulation

Order constraints

Compound effect becomes stronger when dose is increased. Monotone restriction (non-decreasing or non-increasing). Zero effect is meaningful. No parametrical assumptions about dose-response curve shape.

0.8 1.0 1.2 1.4 1.6 1.8 2.0 Dose Response 1 2 3

Martin Otava (I-Biostat, UHasselt) BVS 4 / 20

Model formulation

Basic Model

One-way ANOVA model formulation: Yij = µi + εij i = 0, . . . , K − 1 εij ∼ N(0, σ2) j = 0, 1, 2, . . . , ni = ⇒ necessary to incorporate order constraints. Testing the hypothesis H0 : µ0 = µ1 = µ2 = . . . = µK−1 against ordered alternative (one inequality strict) Hup : µ0 ≤ µ1 ≤ µ2 ≤ . . . ≤ µK−1 Hdn : µ0 ≥ µ1 ≥ µ2 ≥ . . . ≥ µK−1

Martin Otava (I-Biostat, UHasselt) BVS 5 / 20

Bayesian variable selection model

Reformulation of model

New notation (non-decreasing trend): E(Yij) = µi =        µ0, i = 0, µ0 +

i

• ℓ=1

δℓ, i = 1, . . . , K − 1 with priors: µ0 ∼ N(ηµ, σ2

µ),

δi ∼ N(ηδi, σ2

δi)I(0, A),

i = 1, . . . , K − 1. ⇒ δi ≥ 0.

Martin Otava (I-Biostat, UHasselt) BVS 6 / 20

Bayesian variable selection model

Set of all models

0.0 0.5 1.0 1.5

g_0 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_1 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_2 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_3 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_4 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_5 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_6 Dose Response

1 2 3 0.0 0.5 1.0 1.5

g_7 Dose Response

1 2 3

Martin Otava (I-Biostat, UHasselt) BVS 7 / 20

Bayesian variable selection model

Sub-hypotheses

Hup : µ0 ≤ µ1 ≤ µ2 ≤ . . . ≤ µK−1 Model Up: Mean Structure z g0 µ0 = µ1 = µ2 = µ3 (0,0,0) g1 µ0 < µ1 = µ2 = µ3 (1,0,0) g2 µ0 = µ1 < µ2 = µ3 (0,1,0) g3 µ0 < µ1 < µ2 = µ3 (1,1,0) g4 µ0 = µ1 = µ2 < µ3 (0,0,1) g5 µ0 < µ1 = µ2 < µ3 (1,0,1) g6 µ0 = µ1 < µ2 < µ3 (0,1,1) g7 µ0 < µ1 < µ2 < µ3 (1,1,1)

Martin Otava (I-Biostat, UHasselt) BVS 8 / 20

Bayesian variable selection model

Modification to BVS

The distribution of δ is continuous. = ⇒ probability of all models except one equals zero! Instead of only sampling δi we need to select which δi occurs in model. Let be zi indicator of δi occurring in the model. zi =

• 1,

δi is included in the model, 0, δi is not included in the model. ⇒ E(Yij) = µ0 +

i

• ℓ=1

zℓδℓ.

Martin Otava (I-Biostat, UHasselt) BVS 9 / 20

Bayesian variable selection model

BVS model formulation

Basic model: Yij ∼ N(µi, σ2) Modeling of mean: E(Yij) = µi = µ0 +

i

• ℓ=1

zℓδℓ. Priors: µ0 ∼ N(ηµ, σ2

µ),

δi ∼ N(ηδi, σ2

δi)I(0, A),

zi ∼ Bernoulli(πi), Hyper Priors: σ−2 ∼ Γ(10−3, 10−3), ηµ ∼ N(0, 106), σ−2

µ

∼ Γ(10−3, 10−3), ηδi ∼ N(0, 106), σ−2

δi

∼ Γ(10−3, 10−3). πi ∼ U(0, 1).

Martin Otava (I-Biostat, UHasselt) BVS 10 / 20

Results interpretation

Posterior mean of µi

Posterior distribution for all dose-specific means. Use posterior mean of such distribution as our estimation. Connection of Bayesian model averaging. = ⇒ posterior model probabilities are weights. ˆ µBVS =

R

• r=0

wr ˆ µr

Martin Otava (I-Biostat, UHasselt) BVS 11 / 20

Results interpretation

Posterior probability of model

Vector z = (z1, . . . , zK−1) uniquely defines the model. Transformation G(z) = 1 + K−1

i=1 zi 2i−1 =

⇒ unique value for each model. In each MCMC iteration we sample one vector z = (z1, . . . , zK−1). Posterior mean of indicator G(z) = r + 1 translates into posterior probability of the model gr. = ⇒ For posterior probabilities holds: P[G(z) = r + 1|data] = P(gr|data).

Martin Otava (I-Biostat, UHasselt) BVS 12 / 20

Results interpretation

Example: BVS model

Incorporating models with equal means results into less decreasing profile. Posterior means are averages of means of particular models at each MCMC iteration. ˆ µBVS =

R

• r=0

¯ P(gr|data)ˆ µr Connection to model averaging.

• Dose

Weight

1 2 3 20 25 30 35

Order restricted BVS

Martin Otava (I-Biostat, UHasselt) BVS 13 / 20

Results interpretation

Example: Posterior probabilities

Posterior probabilities of particular models. Model g0 represents H0. Model g1 is strongly supported by the data. Connection to model selection.

g_0 g_1 g_2 g_3 g_4 g_5 g_6 g_7 Model Posterior probability 0.0 0.1 0.2 0.3 0.4 0.5 0.6

Martin Otava (I-Biostat, UHasselt) BVS 14 / 20

Hypothesis testing

Hypothesis testing

Depends on: data on hand, prior distributions, set of alternative hypotheses. We use objective priors and consider the set of all possible alternative hypotheses. Use ¯ P(g0|data), estimation of P(H0|data), to reject H0. Questions:

How to select threshold for deciding if H0 is rejected? There is no straightforward control mechanism like Type I error.

Simulation study can give us insight in the properties of BVS.

Martin Otava (I-Biostat, UHasselt) BVS 15 / 20

Hypothesis testing

Simulation study

Under the H0 and under model g7. P(H0|data) < τ used as criterion for rejecting H0 by BVS. PH0(data∗) < τ used as criterion for rejecting H0 by LRT and MCTs. What happens to false rejections and false non-rejections while varying threshold τ? When maintaining approximately same empirical Type I error as MCTs or LRT, BVS seems to achieve similar power. How to select threshold for BVS in practice? = ⇒ future research.

Martin Otava (I-Biostat, UHasselt) BVS 16 / 20

Hypothesis testing

Simulation study - Results

0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0

Threshold Power BVS LRT MCT−W MCT−M

0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0

Threshold Type I error BVS LRT MCT−W MCT−M Martin Otava (I-Biostat, UHasselt) BVS 17 / 20

Discussion

Conclusion

Model uncertainty taken into account! Model selection: ¯ P(gr|data). Estimation of means: ˆ µ = R

r=0 ¯

P(gr|data)ˆ µr. Inference: ¯ P(g0|data). BVS framework address all perspectives simultaneously. According to simulations seems to perform comparably with LRT and MCTs.

Martin Otava (I-Biostat, UHasselt) BVS 18 / 20

Discussion

Future research

How to select threshold for rejecting H0 using P(H0|data)? How to fit BVS models with different types of restrictions (e.g. umbrella profiles)? How do BVS models behave when used for multiplicity adjustment?

Martin Otava (I-Biostat, UHasselt) BVS 19 / 20

Discussion

Thank you for your attention!

Martin Otava (I-Biostat, UHasselt) BVS 20 / 20