Biology Reference
In-Depth Information
> dbn1 =
+ model2network("[245094_at][265768_at|245094_at]")
> xp.mean = mean(x[, "245094_at"])
> xp.sd = sd(x[, "245094_at"])
> dbn1.fit =
+ custom.fit(dbn1,
+ dist = list("245094_at" = list(coef = xp.mean,
+ sd = xp.sd), "265768_at" = lasso.t))
Since we are modeling continuous data, we create a Gaussian Bayesian network. For
the distribution of gene 245094 at , we only need to specify the mean ( xp.mean )
and the standard deviation ( xp.sd ) since the corresponding node has no parents
in dbn1 . For gene 265768 at , we reuse the parameters we estimated for the
lasso.t model.
As expected from the regression coefficient in lasso.t , high expression levels
of gene 245094 at at time t
1 make high expression levels of gene 265768 at
at time t much less likely than lower expression levels of gene 245094 at at time
t
1.
> cpquery(dbn1.fit, event = ('265768_at' > 8),
+
evidence = ('245094_at' > 8))
[1] 0.2448749
> cpquery(dbn1.fit, event = ('265768_at' > 8),
+
evidence = ('245094_at' < 8))
[1] 0.9827778
It is important to note that both event and evidence must specify interval events,
because any single point in
has probability zero. Therefore, any conditional prob-
ability query based on such an event will always return to zero.
For a graphical comparison, we can use cpdist to generate two sets of random
observations under the different conditioning events and compare their densities.
The resulting plot is shown in Fig. 4.3 .
R
> dist.low = cpdist(dbn1.fit, node = "265768_at",
+ evidence = ('245094_at' < 8))
> dist.high = cpdist(dbn1.fit, node = "265768_at",
+ evidence = ('245094_at' > 8))
Performing conditional probability queries using time points that are farther
apart, i.e., s
2and t , is also possible. Consider one more time the expression
of gene 265768 at at time t ; the variables at time t
=
t
2 that are most relevant for
our queries are the parents of gene 245094 at as identified by the LASSO model
having 245094 at as a response variable. To avoid name clashes, the expression
level of gene 245094 at at time t
2 will be called 245094 at1 .
> y = arth12[-(1:2), "245094_at"]
> colnames(x)[12] = "245094_at1"
> lambda = optL1(response = y, penalized = x)$lambda
 
Search WWH ::




Custom Search