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