When learning a hidden Markov model (HMM), sequential observations can often be complemented by real-valued summary response variables generated from the path of hidden states. Such settings arise in numerous domains, including many applications in biology, like motif discovery and genome annotation. In this paper, we present a flexible framework for jointly modeling both latent sequence features and the functional mapping that relates the summary response variables to the hidden state sequence. The algorithm is compatible with a rich set of mapping functions. Results show that the availability of additional continuous response variables can simultaneously improve the annotation of the sequential observations and yield good prediction performance in both synthetic data and real-world datasets.

MOTIVATION: Transcriptional regulation is directly enacted by the interactions between DNA and many proteins, including transcription factors (TFs), nucleosomes and polymerases. A critical step in deciphering transcriptional regulation is to infer, and eventually predict, the precise locations of these interactions, along with their strength and frequency. While recent datasets yield great insight into these interactions, individual data sources often provide only partial information regarding one aspect of the complete interaction landscape. For example, chromatin immunoprecipitation (ChIP) reveals the binding positions of a protein, but only for one protein at a time. In contrast, nucleases like MNase and DNase can be used to reveal binding positions for many different proteins at once, but cannot easily determine the identities of those proteins. Currently, few statistical frameworks jointly model these different data sources to reveal an accurate, holistic view of the in vivo protein-DNA interaction landscape. RESULTS: Here, we develop a novel statistical framework that integrates different sources of experimental information within a thermodynamic model of competitive binding to jointly learn a holistic view of the in vivo protein-DNA interaction landscape. We show that our framework learns an interaction landscape with increased accuracy, explaining multiple sets of data in accordance with thermodynamic principles of competitive DNA binding. The resulting model of genomic occupancy provides a precise mechanistic vantage point from which to explore the role of protein-DNA interactions in transcriptional regulation. AVAILABILITY AND IMPLEMENTATION: The C source code for compete and Python source code for MCMC-based inference are available at http://www.cs.duke.edu/∼amink. CONTACT: [email protected] SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

MOTIVATION: The DNA binding specificity of a transcription factor (TF) is typically represented using a position weight matrix model, which implicitly assumes that individual bases in a TF binding site contribute independently to the binding affinity, an assumption that does not always hold. For this reason, more complex models of binding specificity have been developed. However, these models have their own caveats: they typically have a large number of parameters, which makes them hard to learn and interpret. RESULTS: We propose novel regression-based models of TF-DNA binding specificity, trained using high resolution in vitro data from custom protein-binding microarray (PBM) experiments. Our PBMs are specifically designed to cover a large number of putative DNA binding sites for the TFs of interest (yeast TFs Cbf1 and Tye7, and human TFs c-Myc, Max and Mad2) in their native genomic context. These high-throughput quantitative data are well suited for training complex models that take into account not only independent contributions from individual bases, but also contributions from di- and trinucleotides at various positions within or near the binding sites. To ensure that our models remain interpretable, we use feature selection to identify a small number of sequence features that accurately predict TF-DNA binding specificity. To further illustrate the accuracy of our regression models, we show that even in the case of paralogous TF with highly similar position weight matrices, our new models can distinguish the specificities of individual factors. Thus, our work represents an important step toward better sequence-based models of individual TF-DNA binding specificity. AVAILABILITY: Our code is available at http://genome.duke.edu/labs/gordan/ISMB2013. The PBM data used in this article are available in the Gene Expression Omnibus under accession number GSE47026.

MOTIVATION: To advance understanding of eukaryotic cell division, it is important to observe the process precisely. To this end, researchers monitor changes in dividing cells as they traverse the cell cycle, with the presence or absence of morphological or genetic markers indicating a cell's position in a particular interval of the cell cycle. A wide variety of marker data is available, including information-rich cellular imaging data. However, few formal statistical methods have been developed to use these valuable data sources in estimating how a population of cells progresses through the cell cycle. Furthermore, existing methods are designed to handle only a single binary marker of cell cycle progression at a time. Consequently, they cannot facilitate comparison of experiments involving different sets of markers. RESULTS: Here, we develop a new sampling model to accommodate an arbitrary number of different binary markers that characterize the progression of a population of dividing cells along a branching process. We engineer a strain of Saccharomyces cerevisiae with fluorescently labeled markers of cell cycle progression, and apply our new model to two image datasets we collected from the strain, as well as an independent dataset of different markers. We use our model to estimate the duration of post-cytokinetic attachment between a S.cerevisiae mother and daughter cell. The Java implementation is fast and extensible, and includes a graphical user interface. Our model provides a powerful and flexible cell cycle analysis tool, suitable to any type or combination of binary markers. AVAILABILITY: The software is available from: http://www.cs.duke.edu/~amink/software/cloccs/. CONTACT: [email protected]; [email protected].

MOTIVATION: Recent advances in high-throughput experimental techniques have yielded a large amount of data on protein-protein interactions (PPIs). Since these interactions can be organized into networks, and since separate PPI networks can be constructed for different species, a natural research direction is the comparative analysis of such networks across species in order to detect conserved functional modules. This is the task of network alignment. RESULTS: Most conventional network alignment algorithms adopt a node-then-edge-alignment paradigm: they first identify homologous proteins across networks and then consider interactions among them to construct network alignments. In this study, we propose an alternative direct-edge-alignment paradigm. Specifically, instead of explicit identification of homologous proteins, we directly infer plausibly alignable PPIs across species by comparing conservation of their constituent domain interactions. We apply our approach to detect conserved protein complexes in yeast-fly and yeast-worm PPI networks, and show that our approach outperforms two recent approaches in most alignment performance metrics. AVAILABILITY: Supplementary material and source code can be found at http://www.cs.duke.edu/ approximately amink/.

A principled mechanism for identifying conditional dependencies in time-series data is provided through structure learning of dynamic Bayesian networks (DBNs). An important assumption of DBN structure learning is that the data are generated by a stationary process—an assumption that is not true in many important settings. In this paper, we introduce a new class of graphical models called non-stationary dynamic Bayesian networks, in which the conditional dependence structure of the underlying data-generation process is permitted to change over time. Non-stationary dynamic Bayesian networks represent a new framework for studying problems in which the structure of a network is evolving over time. We define the non-stationary DBN model, present an MCMC sampling algorithm for learning the structure of the model from time-series data under different assumptions, and demonstrate the effectiveness of the algorithm on both simulated and biological data.

Finding functional DNA binding sites of transcription factors (TFs) throughout the genome is a crucial step in understanding transcriptional regulation. Unfortunately, these binding sites are typically short and degenerate, posing a significant statistical challenge: many more matches to known TF motifs occur in the genome than are actually functional. However, information about chromatin structure may help to identify the functional sites. In particular, it has been shown that active regulatory regions are usually depleted of nucleosomes, thereby enabling TFs to bind DNA in those regions. Here, we describe a novel motif discovery algorithm that employs an informative prior over DNA sequence positions based on a discriminative view of nucleosome occupancy. When a Gibbs sampling algorithm is applied to yeast sequence-sets identified by ChIP-chip, the correct motif is found in 52% more cases with our informative prior than with the commonly used uniform prior. This is the first demonstration that nucleosome occupancy information can be used to improve motif discovery. The improvement is dramatic, even though we are using only a statistical model to predict nucleosome occupancy; we expect our results to improve further as high-resolution genome-wide experimental nucleosome occupancy data becomes increasingly available.

MOTIVATION: A key goal in molecular biology is to understand the mechanisms by which a cell regulates the transcription of its genes. One important aspect of this transcriptional regulation is the binding of transcription factors (TFs) to their specific cis-regulatory counterparts on the DNA. TFs recognize and bind their DNA counterparts according to the structure of their DNA-binding domains (e.g. zinc finger, leucine zipper, homeodomain). The structure of these domains can be used as a basis for grouping TFs into classes. Although the structure of DNA-binding domains varies widely across TFs generally, the TFs within a particular class bind to DNA in a similar fashion, suggesting the existence of class-specific features in the DNA sequences bound by each class of TFs. RESULTS: In this paper, we apply a sparse Bayesian learning algorithm to identify a small set of class-specific features in the DNA sequences bound by different classes of TFs; the algorithm simultaneously learns a true multi-class classifier that uses these features to predict the DNA-binding domain of the TF that recognizes a particular set of DNA sequences. We train our algorithm on the six largest classes in TRANSFAC, comprising a total of 587 TFs. We learn a six-class classifier for this training set that achieves 87% leave-one-out cross-validation accuracy. We also identify features within cis-regulatory sequences that are highly specific to each class of TF, which has significant implications for how TF binding sites should be modeled for the purpose of motif discovery.

Determining how information flows along anatomical brain pathways is a fundamental requirement for understanding how animals perceive their environments, learn, and behave. Attempts to reveal such neural information flow have been made using linear computational methods, but neural interactions are known to be nonlinear. Here, we demonstrate that a dynamic Bayesian network (DBN) inference algorithm we originally developed to infer nonlinear transcriptional regulatory networks from gene expression data collected with microarrays is also successful at inferring nonlinear neural information flow networks from electrophysiology data collected with microelectrode arrays. The inferred networks we recover from the songbird auditory pathway are correctly restricted to a subset of known anatomical paths, are consistent with timing of the system, and reveal both the importance of reciprocal feedback in auditory processing and greater information flow to higher-order auditory areas when birds hear natural as opposed to synthetic sounds. A linear method applied to the same data incorrectly produces networks with information flow to non-neural tissue and over paths known not to exist. To our knowledge, this study represents the first biologically validated demonstration of an algorithm to successfully infer neural information flow networks.

MOTIVATION: Duplication of an organism's entire genome is a rare but spectacular event, enabling the rapid emergence of multiple new gene functions. Over time, the parallel linkage of duplicated genes across chromosomes may be disrupted by reciprocal translocations, while the intra-chromosomal order of genes may be shuffled by inversions and transpositions. Some duplicate genes may evolve unrecognizably or be deleted. As a consequence, the only detectable signature of an ancient duplication event in a modern genome may be the presence of various chromosomal segments containing parallel paralogous genes, with each segment appearing exactly twice in the genome. The problem of reconstructing the linkage structure of an ancestral genome before duplication is known as genome halving with unordered chromosomes. RESULTS: In this paper, we derive a new upper bound on the genome halving distance that is tighter than the best known, and a new lower bound that is almost always tighter than the best known. We also define the notion of genome halving diameter, and obtain both upper and lower bounds for it. Our tighter bounds on genome halving distance yield a new algorithm for reconstructing an ancestral duplicated genome. We create a software package GenomeHalving based on this new algorithm and test it on the yeast genome, identifying a sequence of translocations for halving the yeast genome that is shorter than previously conjectured possible.

Recently developed methods for learning sparse classifiers are among the state-of-the-art in supervised learning. These methods learn classifiers that incorporate weighted sums of basis functions with sparsity-promoting priors encouraging the weight estimates to be either significantly large or exactly zero. From a learning-theoretic perspective, these methods control the capacity of the learned classifier by minimizing the number of basis functions used, resulting in better generalization. This paper presents three contributions related to learning sparse classifiers. First, we introduce a true multiclass formulation based on multinomial logistic regression. Second, by combining a bound optimization approach with a component-wise update procedure, we derive fast exact algorithms for learning sparse multiclass classifiers that scale favorably in both the number of training samples and the feature dimensionality, making them applicable even to large data sets in high-dimensional feature spaces. To the best of our knowledge, these are the first algorithms to perform exact multinomial logistic regression with a sparsity-promoting prior. Third, we show how nontrivial generalization bounds can be derived for our classifier in the binary case. Experimental results on standard benchmark data sets attest to the accuracy, sparsity, and efficiency of the proposed methods.

A graph-based prior is proposed for parametric semi-supervised classi- fication. The prior utilizes both labelled and unlabelled data; it also in- tegrates features from multiple views of a given sample (e.g., multiple sensors), thus implementing a Bayesian form of co-training. An EM algorithm for training the classifier automatically adjusts the tradeoff be- tween the contributions of: (a) the labelled data; (b) the unlabelled data; and (c) the co-training information. Active label query selection is per- formed using a mutual information based criterion that explicitly uses the unlabelled data and the co-training information. Encouraging results are presented on public benchmarks and on measured data from single and multiple sensors. 1 Introduction In many pattern classification problems, the acquisition of labelled training data is costly and/or time consuming, whereas unlabelled samples can be obtained easily. Semi- supervised algorithms that learn from both labelled and unlabelled samples have been the focus of much research in the last few years; a comprehensive review up to 2001 can be found in [13], while more recent references include [1, 2, 6, 7, 1618]. Most recent semi-supervised learning algorithms work by formulating the assumption that "nearby" points, and points in the same structure (e.g., cluster), should have similar labels [6, 7, 16]. This can be seen as a form of regularization, pushing the class boundaries toward regions of low data density. This regularization is often implemented by associating the vertices of a graph to all the (labelled and unlabelled) samples, and then formulating the problem on the vertices of the graph [6, 1618]. While current graph-based algorithms are inherently transductive -- i.e., they cannot be used directly to classify samples not present when training -- our classifier is paramet- ric and the learned classifier can be used directly on new samples. Furthermore, our al- gorithm is trained discriminatively by maximizing a concave objective function; thus we avoid thorny local maxima issues that plague many earlier methods. Unlike existing methods, our algorithm automatically learns the relative importance of the labelled and unlabelled data. When multiple views of the same sample are provided (e.g. features from different sensors), we develop a new Bayesian form of co-training [4]. In addition, we also show how to exploit the unlabelled data and the redundant views of the sample (from co-training) in order to improve active label query selection [15]. The paper is organized as follows. Sec. 2 briefly reviews multinomial logistic regression. Sec. 3 describes the priors for semi-supervised learning and co-training. The EM algorithm derived to learn the classifiers is presented in Sec. 4. Active label selection is discussed in Sec. 5. Experimental results are shown in Sec. 6, followed by conclusions in Sec. 7. 2 Multinomial Logistic Regression In an m-class supervised learning problem, one is given a labelled training set DL = {(x d 1, y ), . . . , (x )} 1 L, yL , where xi R is a feature vector and yi the corresponding (1) (m) class label. In "1-of-m" encoding, y = [y , . . . , y ] i is a binary vector, such that i i (c) (j) y = 1 and y = 0, for j = c, indicates that sample i belongs to class c. In multinomial i i logistic regression [5], the posterior class probabilities are modelled as log P (y(c) = 1|x) = xT w(c) - log m exp(xT w(k)), for c = 1, . . . , m, (1) k=1 where w(c) d R is the class-c weight vector. Notice that since m P (y(c)= 1|x) = 1, c=1 one of the weight vectors is redundant; we arbitrarily choose to set w(m) = 0, and consider the (d (m-1))-dimensional vector w = [(w(1))T , ..., (w(m-1))T ]T . Estimation of w may be achieved by maximizing the log-likelihood (with Y {y , ..., y } 1 L ) [5] (c) (w) log P (Y|w) = L m y xT w(c) - log m exp(xT w(j)) . (2) i=1 c=1 i i j=1 i In the presence of a prior p(w), we seek a maximum a posteriori (MAP) estimate, w = arg max { (w) + log p(w)} w . Actually, if the training data is separable, (w) is unbounded, and a prior is crucial. Although we focus on linear classifiers, we may see the d-dimensional feature vectors x as having resulted from some deterministic, maybe nonlinear, transformation of an input raw feature vector r; e.g., in a kernel classifier, xi = [1, K(ri, r1), ..., K(ri, rL)] (d = L + 1). 3 Graph-Based Data-Dependent Priors 3.1 Graph Laplacians and Regularization for Semi-Supervised Learning Consider a scalar function f = [f1, ..., f|V |]T , defined on the set V = {1, 2, ..., |V |} of vertices of an undirected graph (V, E). Each edge of the graph, joining vertices i and j, is given a weight kij = kji 0, and we collect all the weights in a |V | |V | matrix K. A natural way to measure how much f varies across the graph is by the quantity kij(fi - fj)2 = 2 f T f , (3) i j where = diag{ k k j 1j , ..., j |V |j } - K is the so-called graph Laplacian [2]. Notice that kij 0 (for all i, j) guarantees that is positive semi-definite and also that has (at least) one null eigenvalue (1T1 = 0, where 1 has all elements equal to one). In semi-supervised learning, in addition to DL, we are given U unlabelled samples DU = {xL+1, . . . , xL+U }. To use (3) for semi-supervised learning, the usual choice is to assign one vertex of the graph to each sample in X = [x1, . . . , xL+U ]T (thus |V | = L + U ), and to let kij represent some (non-negative) measure of "similarity" between xi and xj. A Gaussian random field (GRF) is defined on the vertices of V (with inverse variance ) p(f ) exp{- f T f /2}, in which configurations that vary more (according to (3)) are less probable. Most graph- based approaches estimate the values of f , given the labels, using p(f ) (or some modifica- tion thereof) as a prior. Accordingly, they work in a strictly transductive manner. 3.2 Non-Transductive Semi-Supervised Learning We first consider two-class problems (m = 2, thus w d R ). In contrast to previous uses of graph-based priors, we define f as the real function f (defined over the entire observation space) evaluated at the graph nodes. Specifically, f is defined as a linear function of x, and at the graph node i, fi f (xi) = wT xi. Then, f = [f1, ..., f|V |]T = Xw, and p(f ) induces a Gaussian prior on w, with precision matrix A = XT X, p(w) exp{-(/2) wT XT Xw} = exp{-(/2) wT Aw}. (4) Notice that since is singular, A may also be singular, and the corresponding prior may therefore be improper. This is no problem for MAP estimation of w because (as is well known) the normalization factor of the prior plays no role in this estimate. If we include extra regularization, by adding a non-negative diagonal matrix to A, the prior becomes p(w) exp -(1/2) wT (0A + ) w , (5) where we may choose = diag{1, ..., d}, = 1I, or even = 0. For m > 2, we define (m-1) identical independent priors, one for each w(c), c = 1, ..., m. The joint prior on w = [(w(1))T , ..., (w(m-1))T ]T is then m-1 1 (c) 1 p(w|) exp{- (w(c))T A + (c) w(c)} = exp{- wT ()w}, (6) 2 0 2 c=1 (c) (c) (c) where is a vector containing all the parameters, (c) = diag{ , ..., }, and i 1 d (1) (m-1) () = diag{ , ..., } A + block-diag{(1), ..., (m-1)}. (7) 0 0 Finally, since all the 's are inverses of variances, the conjugate priors are Gamma [3]: (c) (c) (c) (c) p( | | | | 0 0, 0) = Ga(0 0, 0), and p(i 1, 1) = Ga(i 1, 1), for c = 1, ..., m - 1 and i = 1, ..., d. Usually, 0, 0, 1, and 1 are given small values indicating diffuse priors. In the zero limit, we obtain scale-invariant (improper) Jeffreys hyper-priors. Summarizing, our model for semi-supervised learning includes the log-likelihood (2), a prior (6), and Gamma hyper-priors. In Section 4, we present a simple and computationally efficient expectation-maximization (EM) algorithm for obtaining the MAP estimate of w. 3.3 Exploiting Features from Multiple Sensors: The Co-Training Prior In some applications several sensors are available, each providing a different set of features. For simplicity, we assume two sensors s {1, 2}, but everything discussed here is easily (s) extended to any number of sensors. Denote the features from sensor s, for sample i, as x , i and Ss as the set of sample indices for which we have features from sensor s (S1 S2 = {1, ..., L + U }). Let O = S1 S2 be the indices for which both sensors are available, and OU = O {L + 1, ..., L + U } the unlabelled subset of O. By using the samples in S1 and S2 as two independent training sets, we may obtain two sep- arate classifiers (denoted w1 and w2). However, we can coordinate the information from both sensors by using an idea known as co-training [4]: on the OU samples, classifiers w1 and w2 should agree as much as possible. Notice that, in a logistic regression framework, the disagreement between the two classifiers on the OU samples can be measured by (1) (2) [(w1)T x - (w2)T x ]2 = T C , (8) iOU i i where = [(w1)T (w2)T ]T and C = [(x1)T (-x2)T ]T [(x1)T (-x2)T ]. This iOU i i i i suggests the "co-training prior" (where co is an inverse variance): p(w1, w2) = p() exp -(co/2) TC . (9) This Gaussian prior can be combined with two smoothness Gaussian priors on w1 and w2 (obtained as described in Section 3.2); this leads to a prior which is still Gaussian, p(w1, w2) = p() exp -(1/2) T coC + block-diag{1, 2} , (10) where 1 and 2 are the two graph-based precision matrices (see (7)) for w1 and w2. We can again adopt a Gamma hyper-prior for co. Under this prior, and with a logistic regression likelihood as above, estimates of w1 and w2 can easily be found using minor modifications to the EM algorithm described in Section 4. Computationally, this is only slightly more expensive than separately training the two classifiers. 4 Learning Via EM To find the MAP estimate w, we use the EM algorithm, with as missing data, which is equivalent to integrating out from the full posterior before maximization [8]. For simplicity, we will only describe the single sensor case (no co-training). E-step: We compute the expected value of the complete log-posterior, given Y and the current parameter estimate w: Q(w|w) E[log p(w, |Y)|w]. Since log p(w, |Y) = log p(Y|w) - (1/2)wT ()w + K, (11) (where K collects all terms independent of w) is linear w.r.t. all the parameters (see (6) and (7)), we just have to plug their conditional expectations into (11): Q(w|w) = log p(Y|w) - (1/2)wT E[()|w] w = (w) - (1/2)wT (w) w. (12) We consider several different choices for the structure of the matrix. The necessary expectations have well-known closed forms, due to the use of conjugate Gamma hyper- (c) priors [3]. For example, if the are m - 1 free non-negative parameters, we have 0 (c) (c) E[ |w] = (2 0 0 0 + d) [2 0 + (w(c))T Aw(c)]-1. (c) for c = 1, ..., m - 1. For = 0 0, we still have a simple closed-form expres- (c) sion for E[0|w], and the same is true for the parameters, for i > 0. Finally, i (w) E[()|w] results from replacing the 's in (7) by the corresponding conditional expectations. M-step: Given matrix (w), the M-step reduces to a logistic regression problem with a quadratic regularizer, i.e., maximizing (12). To this end, we adopt the bound optimization approach (see details in [5, 11]). Let B be a positive definite matrix such that -B bounds below (in the matrix sense) the Hessian of (w), which is negative definite, and g(w) is the gradient of (w). Then, we have the following lower bound on Q(w|w): Q(w|w) l(w) + (w - w)T g(w) - [(w - w)T B(w - w) + wT (w)w]/2. - The maximizer of this lower bound, wnew = (B + (w)) 1 (Bw + g(w)), is guaranteed to increase the Q-function, Q(wnew|w) Q(w|w), and we thus obtain a monotonic gen- eralized EM algorithm [5, 11]. This (maybe costly) matrix inversion can be avoided by a sequential approach where we only maximize w.r.t. one element of w at a time, preserving the monotonicity of the procedure. The sequential algorithm visits one particular element of w, say wu, and updates its estimate by maximizing the bound derived above, while keeping all other variables fixed at their previous values. This leads to - wnew = w ] [(B + (w)) 1 , u u + [gu(w) - ((w)w) (13) u uu] and wnew = w v v , for v = u. The total time required by a full sweep for all u = 1, ..., d is O(md(L + d)); this may be much better than the O((dm)3) of the matrix inversion. 5 Active Label Selection If we are allowed to obtain the label for one of the unlabelled samples, the following ques- tion arises: which sample, if labelled, would provide the most information? Consider the MAP estimate w provided by EM. Our approach uses a Laplace approxima- tion of the posterior p(w|Y) N (w|w, H-1), where H is the posterior precision matrix, i.e., the Hessian of minus the log-posterior H = 2(- log p(w|Y)). This approximation is known to be accurate for logistic regression under a Gaussian prior [14]. By treating (w) (the expectation of ()) as deterministic, we obtain an evidence-type approximation [14] H = 2[- log(p(Y|w)p(w|(w)))] = (w) + L (diag{p } - p pT ) x , i=1 i i i ixT i where pi is the (m - 1)-dimensional vector computed from (1), the c-th element of which indicates the probability that sample xi belongs to class c. Now let x DU be an unlabelled sample and y its label. Assume that the MAP esti- mate w remains unchanged after including y. In Sec. 7 we will discuss the merits and shortcomings of this assumption, which is only strictly valid when L . Accepting it implies that after labeling x, and regardless of y, the posterior precision changes to H = H + (diag{p} - ppT ) xxT . (14) Since the entropy of a Gaussian with precision H is (-1/2) log |H| (up to an additive constant), the mutual information (MI) between y and w (i.e., the expected decrease in entropy of w when y is observed) is I(w; y) = (1/2) log {|H |/|H|}. Our criterion is then: the best sample to label is the one that maximizes I(w; y). Further insight into I(w; y) can be obtained in the binary case (where p is a scalar); here, the matrix identity |H + p(1 - p)xxT | = |H|(1 + p(1 - p)xT H-1x) yields I(w; y) = (1/2) log(1 + p(1 - p)xT H-1x). (15) This MI is larger when p 0.5, i.e., for samples with uncertain classifications. On the other hand, with p fixed, I(w; y) grows with xT H-1x, i.e., it is large for samples with high variance of the corresponding class probability estimate. Summarizing, (15) favors samples with uncertain class labels and high uncertainty in the class probability estimate. 6 Experimental Results We begin by presenting two-dimensional synthetic examples to visually illustrate our semi- supervised classifier. Fig. 1 shows the utility of using unlabelled data to improve the deci- Figure 1: Synthetic two-dimensional examples. (a) Comparison of the supervised logistic linear classifier (boundary shown as dashed line) learned only from the labelled data (shown in color) with the proposed semi-supervised classifier (boundary shown as solid line) which also uses the unlabelled samples (shown as dots). (b) A RBF kernel classifier obtained by our algorithm, using two labelled samples (shaded circles) and many unlabelled samples. Figure 2: (a)-(c) Accuracy (on UCI datasets) of the proposed method, the supervised SVM, and the other semi-supervised classifiers mentioned in the text; a subset of samples is la- belled and the others are treated as unlabelled samples. In (d), a separate holdout set is used to evaluate the accuracy of our method versus the amount of labelled and unlabelled data. sion boundary in linear and non-linear (kernel) classifiers (see figure caption for details). Next we show results with linear classifiers on three UCI benchmark datasets. Results with nonlinear kernels are similar, and therefore omitted to save space. We compare our method against state-of-the-art semi-supervised classifiers: the GRF method of [18], the SGT method of [10], and the transductive SVM (TSVM) of [9]. For reference, we also present results for a standard SVM. To avoid unduly helping our method, we always use a k=5 nearest neighbors graph, though our algorithm is not very sensitive to k. To avoid disadvantaging other methods that do depend on such parameters, we use their best settings. Since these adjustments cannot be made in practice, the difference between our algorithm and the others is under-represented. Each point on the plots in Fig. 2(a)-(c) is an average of 20 trials: we randomly select 20 labelled sets which are used by every method. All remaining samples are used as unlabelled by the semi-supervised algorithms. Figs. 2(a)-(c) are transductive, in the sense that the unlabelled and test data are the same. Our logistic GRF is non-transductive: after being trained, it may be applied to classify new data without re-training. In Fig. 2(d) we present non-transductive results for the Ionosphere data. Training took place using labelled and unlabelled data, and testing was performed on 200 new unseen samples. The results suggest that semi-supervised classifiers are most relevant when the labelled set is small relative to the unlabelled set (as is often the case). Our final set of results address co-training (Sec. 3.3) and active learning (Sec. 5), applied to airborne sensing data for the detection of surface and subsurface land mines. Two sensors were used: (1) a 70-band hyper-spectral electro-optic (EOIR) sensor; (2) an X-band syn- thetic aperture radar (SAR). A simple (energy) "prescreener" detected potential targets; for each of these, two feature vectors were extracted, of sizes 420 and 9, for the EOIR and SAR sensors, respectively. 123 samples have features from the EOIR sensor alone, 398 from the Figure 3: (a) Land mine detection ROC curves of classifiers designed using only hyper- spectral (EOIR) features, only SAR features, and both. (b) Number of landmines detected during the active querying process (dotted lines), for active training and random selection (for the latter the bars reflect one standard deviation about the mean). ROC curves (solid) are for the learned classifier as applied to the remaining samples. SAR sensor alone, and 316 from both. This data will be made available upon request. We first consider supervised and semi-supervised classification. For the purely supervised case, a sparseness prior is used (as in [14]). In both cases a linear classifier is employed. For the data for which only one sensor is available, 20% of it is labelled (selected randomly). For the data for which both sensors are available, 80% is labelled (again selected randomly). The results presented in Fig. 3(a) show that, in general, the semi-supervised classifiers outperform the corresponding supervised ones, and the classifier learned from both sensors is markedly superior to classifiers learned from either sensor alone. In a second illustration, we use the active-learning algorithm (Sec. 5) to only acquire the 100 most informative labels. For comparison, we also show average results over 100 in- dependent realizations for random label query selection (error bars indicate one standard deviation). The results in Fig. 3(b) are plotted in two stages: first, mines and clutter are se- lected during the labeling process (dashed curves); then, the 100 labelled examples are used to build the final semi-supervised classifier, for which the ROC curve is obtained using the remaining unlabelled data (solid curves). Interestingly, the active-learning algorithm finds almost half of the mines while querying for labels. Due to physical limitations of the sen- sors, the rate at which mines are detected drops precipitously after approximately 90 mines are detected -- i.e., the remaining mines are poorly matched to the sensor physics.

MOTIVATION: Network inference algorithms are powerful computational tools for identifying putative causal interactions among variables from observational data. Bayesian network inference algorithms hold particular promise in that they can capture linear, non-linear, combinatorial, stochastic and other types of relationships among variables across multiple levels of biological organization. However, challenges remain when applying these algorithms to limited quantities of experimental data collected from biological systems. Here, we use a simulation approach to make advances in our dynamic Bayesian network (DBN) inference algorithm, especially in the context of limited quantities of biological data. RESULTS: We test a range of scoring metrics and search heuristics to find an effective algorithm configuration for evaluating our methodological advances. We also identify sampling intervals and levels of data discretization that allow the best recovery of the simulated networks. We develop a novel influence score for DBNs that attempts to estimate both the sign (activation or repression) and relative magnitude of interactions among variables. When faced with limited quantities of observational data, combining our influence score with moderate data interpolation reduces a significant portion of false positive interactions in the recovered networks. Together, our advances allow DBN inference algorithms to be more effective in recovering biological networks from experimentally collected data. AVAILABILITY: Source code and simulated data are available upon request. SUPPLEMENTARY INFORMATION: http://www.jarvislab.net/Bioinformatics/BNAdvances/

This paper adopts a Bayesian approach to simultaneously learn both an optimal nonlinear classifier and a subset of predictor variables (or features) that are most relevant to the classification task. The approach uses heavy-tailed priors to promote sparsity in the utilization of both basis functions and features; these priors act as regularizers for the likelihood function that rewards good classification on the training data. We derive an expectation-maximization (EM) algorithm to efficiently compute a maximum a posteriori (MAP) point estimate of the various parameters. The algorithm is an extension of recent state-of-the-art sparse Bayesian classifiers, which in turn can be seen as Bayesian counterparts of support vector machines. Experimental comparisons using kernel classifiers demonstrate both parsimonious feature selection and excellent classification accuracy on a range of synthetic and benchmark data sets.

Recent research has demonstrated quite convincingly that accurate cancer diagnosis can be achieved by constructing classifiers that are designed to compare the gene expression profile of a tissue of unknown cancer status to a database of stored expression profiles from tissues of known cancer status. This paper introduces the JCFO, a novel algorithm that uses a sparse Bayesian approach to jointly identify both the optimal nonlinear classifier for diagnosis and the optimal set of genes on which to base that diagnosis. We show that the diagnostic classification accuracy of the proposed algorithm is superior to a number of current state-of-the-art methods in a full leave-one-out cross-validation study of two widely used benchmark datasets. In addition to its superior classification accuracy, the algorithm is designed to automatically identify a small subset of genes (typically around twenty in our experiments) that are capable of providing complete discriminatory information for diagnosis. Focusing attention on a small subset of genes is not only useful because it produces a classifier with good generalization capacity, but also because this set of genes may provide insights into the mechanisms responsible for the disease itself. A number of the genes identified by the JCFO in our experiments are already in use as clinical markers for cancer diagnosis; some of the remaining genes may be excellent candidates for further clinical investigation. If it is possible to identify a small set of genes that is indeed capable of providing complete discrimination, inexpensive diagnostic assays might be widely deployable in clinical settings.

MOTIVATION: Although many network inference algorithms have been presented in the bioinformatics literature, no suitable approach has been formulated for evaluating their effectiveness at recovering models of complex biological systems from limited data. To overcome this limitation, we propose an approach to evaluate network inference algorithms according to their ability to recover a complex functional network from biologically reasonable simulated data. RESULTS: We designed a simulator to generate data representing a complex biological system at multiple levels of organization: behaviour, neural anatomy, brain electrophysiology, and gene expression of songbirds. About 90% of the simulated variables are unregulated by other variables in the system and are included simply as distracters. We sampled the simulated data at intervals as one would sample from a biological system in practice, and then used the sampled data to evaluate the effectiveness of an algorithm we developed for functional network inference. We found that our algorithm is highly effective at recovering the functional network structure of the simulated system-including the irrelevance of unregulated variables-from sampled data alone. To assess the reproducibility of these results, we tested our inference algorithm on 50 separately simulated sets of data and it consistently recovered almost perfectly the complex functional network structure underlying the simulated data. To our knowledge, this is the first approach for evaluating the effectiveness of functional network inference algorithms at recovering models from limited data. Our simulation approach also enables researchers a priori to design experiments and data-collection protocols that are amenable to functional network inference.