DAG construction

Introduction

The posterior probabilities computed by the findr function can be used as a graph structure prior in Bayesian network learning using the model introduced in the paper High-dimensional Bayesian network inference from systems genetics data using genetic node ordering. This involves converting the pairwise posterior probabilities to a directed acyclic graph (DAG).

BioFindr implements the original “greedy edges” algorithm where edges are added one-by-one in decreasing order of probability, and only if they do not create a cycle in the graph, using an incremental cycle detection algorithm. Two additional algorithms from the paper Maximal acyclic subgraph optimization for gene regulatory networks are also implemented: the heuristic sort algorithm where vertices are sorted by their ratio of out-degree to in-degree, and edges are added only if their source vertex precedes their target vertex in the sorted list, and the greedy insertion algorithm where vertices are iteratively inserted in the position in the current ordering that yields the maximum possible gain of edge weights, where the gain is counted as the difference between the sum of new edge weights included and the sum of old edge weights lost, and edges are counted only if their source vertex precedes their target vertex in the ordering.

Set up the environment

using DrWatson
using DataFrames
using Arrow
using GLMakie
using GraphMakie: graphplot

using BioFindr

Load data

For illustration we use the GEUVADIS microRNA data:

dm = DataFrame(Arrow.Table(datadir("exp_pro","findr-data-geuvadis", "dm.arrow")));
dgm = DataFrame(Arrow.Table(datadir("exp_pro","findr-data-geuvadis", "dgm.arrow")));
dpm = DataFrame(SNP_ID = names(dgm), GENE_ID=names(dm)[1:ncol(dgm)]);

Run BioFindr

We perform causal inference to compute posterior probabilities from all microRNAs with an eQTL to the total set of microRNAs:

dP = findr(dm, dgm, dpm; FDR=0.25)
55×4 DataFrame
30 rows omitted
Row Source Target Probability qvalue
String String Float64 Float64
1 hsa-miR-550a-2-5p hsa-miR-550a-3-5p 1.0 0.0
2 hsa-miR-550a-2-3p hsa-miR-550a-3-5p 1.0 1.26598e-10
3 hsa-miR-3130-1-5p hsa-miR-3130-2-5p 1.0 5.15569e-10
4 hsa-miR-550a-2-3p hsa-miR-550a-2-5p 1.0 8.94978e-10
5 hsa-miR-550a-2-5p hsa-miR-550a-2-3p 1.0 1.4715e-9
6 hsa-miR-335-5p hsa-miR-335-3p 0.999995 8.15308e-7
7 hsa-miR-335-3p hsa-miR-335-5p 0.999921 1.19308e-5
8 hsa-miR-574-5p hsa-miR-574-3p 0.999919 2.05171e-5
9 hsa-miR-3130-2-5p hsa-miR-3130-1-5p 0.987371 0.00142145
10 hsa-miR-1908-3p hsa-miR-1908-5p 0.938239 0.00745537
11 hsa-miR-3667-5p hsa-miR-3667-3p 0.898732 0.0159838
12 hsa-miR-3176-3p hsa-miR-19b-1-3p 0.893081 0.0235617
13 hsa-miR-3667-3p hsa-miR-3667-5p 0.88392 0.0306785
⋮ ⋮ ⋮ ⋮ ⋮
44 hsa-miR-641-5p hsa-miR-103a-2-5p 0.636497 0.214726
45 hsa-miR-3176-3p hsa-miR-130b-3p 0.634758 0.218071
46 hsa-miR-1304-3p hsa-miR-619-5p 0.633773 0.221292
47 hsa-miR-1307-5p hsa-miR-589-5p 0.63287 0.224395
48 hsa-miR-3176-3p hsa-miR-26a-1-5p 0.627556 0.227479
49 hsa-miR-3176-3p hsa-miR-744-3p 0.627122 0.230447
50 hsa-miR-641-5p hsa-miR-103a-1-5p 0.624291 0.233352
51 hsa-miR-193b-3p hsa-miR-615-3p 0.624061 0.236148
52 hsa-miR-1304-3p hsa-miR-16-2-3p 0.615244 0.239005
53 hsa-miR-3176-3p hsa-miR-143-3p 0.601226 0.24202
54 hsa-miR-641-5p hsa-miR-26a-2-3p 0.599607 0.244953
55 hsa-miR-3176-3p hsa-miR-4524a-5p 0.594893 0.247865

To construct a DAG from the DataFrame dP with posterior probabilities using the default “greedy edges” algorithm, run:

G, name2idx = dagfindr!(dP);

The ! in the dagfindr! function name indicates that the function modifies its input argument. We can see that dP indeed contains some new columns:

dP
55×7 DataFrame
30 rows omitted
Row Source Target Probability qvalue Source_idx Target_idx inDAG_greedy_edges
String String Float64 Float64 Int64 Int64 Bool
1 hsa-miR-550a-2-5p hsa-miR-550a-3-5p 1.0 0.0 1 18 true
2 hsa-miR-550a-2-3p hsa-miR-550a-3-5p 1.0 1.26598e-10 2 18 true
3 hsa-miR-3130-1-5p hsa-miR-3130-2-5p 1.0 5.15569e-10 3 7 true
4 hsa-miR-550a-2-3p hsa-miR-550a-2-5p 1.0 8.94978e-10 2 1 true
5 hsa-miR-550a-2-5p hsa-miR-550a-2-3p 1.0 1.4715e-9 1 2 false
6 hsa-miR-335-5p hsa-miR-335-3p 0.999995 8.15308e-7 4 5 true
7 hsa-miR-335-3p hsa-miR-335-5p 0.999921 1.19308e-5 5 4 false
8 hsa-miR-574-5p hsa-miR-574-3p 0.999919 2.05171e-5 6 19 true
9 hsa-miR-3130-2-5p hsa-miR-3130-1-5p 0.987371 0.00142145 7 3 false
10 hsa-miR-1908-3p hsa-miR-1908-5p 0.938239 0.00745537 8 20 true
11 hsa-miR-3667-5p hsa-miR-3667-3p 0.898732 0.0159838 9 11 true
12 hsa-miR-3176-3p hsa-miR-19b-1-3p 0.893081 0.0235617 10 21 true
13 hsa-miR-3667-3p hsa-miR-3667-5p 0.88392 0.0306785 11 9 false
⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮
44 hsa-miR-641-5p hsa-miR-103a-2-5p 0.636497 0.214726 17 49 true
45 hsa-miR-3176-3p hsa-miR-130b-3p 0.634758 0.218071 10 50 true
46 hsa-miR-1304-3p hsa-miR-619-5p 0.633773 0.221292 13 51 true
47 hsa-miR-1307-5p hsa-miR-589-5p 0.63287 0.224395 16 52 true
48 hsa-miR-3176-3p hsa-miR-26a-1-5p 0.627556 0.227479 10 53 true
49 hsa-miR-3176-3p hsa-miR-744-3p 0.627122 0.230447 10 54 true
50 hsa-miR-641-5p hsa-miR-103a-1-5p 0.624291 0.233352 17 55 true
51 hsa-miR-193b-3p hsa-miR-615-3p 0.624061 0.236148 12 56 true
52 hsa-miR-1304-3p hsa-miR-16-2-3p 0.615244 0.239005 13 57 true
53 hsa-miR-3176-3p hsa-miR-143-3p 0.601226 0.24202 10 58 true
54 hsa-miR-641-5p hsa-miR-26a-2-3p 0.599607 0.244953 17 59 true
55 hsa-miR-3176-3p hsa-miR-4524a-5p 0.594893 0.247865 10 60 true

The Source_idx and Target_idx are numerical IDs for the Source and Target nodes, respectively, and the inDAG_greedy_edges indicates whether the edge represented by a row of dP in included in the output DAG G. The mapping from node names to IDs is also returned as a dictionary object name2idx.

name2idx
Dict{Symbol, Int64} with 60 entries:
  Symbol("hsa-miR-1307-5p")   => 16
  Symbol("hsa-miR-769-5p")    => 42
  Symbol("hsa-miR-3130-2-5p") => 7
  Symbol("hsa-let-7b-5p")     => 14
  Symbol("hsa-miR-320b-2-3p") => 39
  Symbol("hsa-miR-143-3p")    => 58
  Symbol("hsa-miR-548ac-5p")  => 23
  Symbol("hsa-miR-625-5p")    => 30
  Symbol("hsa-miR-3176-3p")   => 10
  Symbol("hsa-miR-26a-1-5p")  => 53
  Symbol("hsa-miR-3667-5p")   => 9
  Symbol("hsa-miR-574-3p")    => 19
  Symbol("hsa-miR-550a-2-3p") => 2
  Symbol("hsa-miR-877-5p")    => 45
  Symbol("hsa-miR-26a-2-3p")  => 59
  Symbol("hsa-miR-582-3p")    => 43
  Symbol("hsa-miR-103a-2-5p") => 49
  Symbol("hsa-miR-1307-3p")   => 15
  Symbol("hsa-miR-3130-1-5p") => 3
  ⋮                           => ⋮

The output G is a directed graph object from the Graphs package:

G
{60, 51} directed simple Int64 graph

This is a fairly simple datastructure, which only supports numerical node IDs, hence the need to create the name2idx map. One useful thing one can do with a Graph object is to draw it:

graphplot(G)

For more details, see the documentation.

To run the dagfindr! function with the other DAG construction algorithms mentioned in the Introduction,

G, name2idx = dagfindr!(dP; method="heuristic sort");
G, name2idx = dagfindr!(dP; method="greedy insertion");

The results of these dagfindr! calls are added to dP allowing easy comparison of the methods:

dP
55×9 DataFrame
30 rows omitted
Row Source Target Probability qvalue Source_idx Target_idx inDAG_greedy_edges inDAG_heuristic_sort inDAG_greedy_insertion
String String Float64 Float64 Int64 Int64 Bool Bool Bool
1 hsa-miR-550a-2-5p hsa-miR-550a-3-5p 1.0 0.0 1 18 true true true
2 hsa-miR-550a-2-3p hsa-miR-550a-3-5p 1.0 1.26598e-10 2 18 true true true
3 hsa-miR-3130-1-5p hsa-miR-3130-2-5p 1.0 5.15569e-10 3 7 true true true
4 hsa-miR-550a-2-3p hsa-miR-550a-2-5p 1.0 8.94978e-10 2 1 true true true
5 hsa-miR-550a-2-5p hsa-miR-550a-2-3p 1.0 1.4715e-9 1 2 false false false
6 hsa-miR-335-5p hsa-miR-335-3p 0.999995 8.15308e-7 4 5 true true true
7 hsa-miR-335-3p hsa-miR-335-5p 0.999921 1.19308e-5 5 4 false false false
8 hsa-miR-574-5p hsa-miR-574-3p 0.999919 2.05171e-5 6 19 true true true
9 hsa-miR-3130-2-5p hsa-miR-3130-1-5p 0.987371 0.00142145 7 3 false false false
10 hsa-miR-1908-3p hsa-miR-1908-5p 0.938239 0.00745537 8 20 true true true
11 hsa-miR-3667-5p hsa-miR-3667-3p 0.898732 0.0159838 9 11 true true true
12 hsa-miR-3176-3p hsa-miR-19b-1-3p 0.893081 0.0235617 10 21 true true true
13 hsa-miR-3667-3p hsa-miR-3667-5p 0.88392 0.0306785 11 9 false false false
⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮ ⋮
44 hsa-miR-641-5p hsa-miR-103a-2-5p 0.636497 0.214726 17 49 true true true
45 hsa-miR-3176-3p hsa-miR-130b-3p 0.634758 0.218071 10 50 true true true
46 hsa-miR-1304-3p hsa-miR-619-5p 0.633773 0.221292 13 51 true true true
47 hsa-miR-1307-5p hsa-miR-589-5p 0.63287 0.224395 16 52 true true true
48 hsa-miR-3176-3p hsa-miR-26a-1-5p 0.627556 0.227479 10 53 true true true
49 hsa-miR-3176-3p hsa-miR-744-3p 0.627122 0.230447 10 54 true true true
50 hsa-miR-641-5p hsa-miR-103a-1-5p 0.624291 0.233352 17 55 true true true
51 hsa-miR-193b-3p hsa-miR-615-3p 0.624061 0.236148 12 56 true true true
52 hsa-miR-1304-3p hsa-miR-16-2-3p 0.615244 0.239005 13 57 true true true
53 hsa-miR-3176-3p hsa-miR-143-3p 0.601226 0.24202 10 58 true true true
54 hsa-miR-641-5p hsa-miR-26a-2-3p 0.599607 0.244953 17 59 true true true
55 hsa-miR-3176-3p hsa-miR-4524a-5p 0.594893 0.247865 10 60 true true true