<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Directed Graphical Models
DGM is a class implementing Directed Graphical Models.
Loading DGMs
You can load bayesian networks in .bif format. For example, here we load a simple 'Earthquake' network.
Some networks can be downloaded here
Step1:  DGM is a subclass of networkx.DiGraph, so you can use any networkx functions on it.
Step2:  Also, some DGM-specific queries about the graph are implemented
Step3:  Checks for I-equivalence and I-maps
Step4:  Find which nodes are reachable in the sense of d-separation
Step5:  Factors
Currently, only TableFactor class is implemented. TableFactor represens an unnormalized discrete CPD.
Step6:  Operations
All factor operations are either vectorized or written in C++, so they should be quite fast.
Step7:  Factor Multiplication
Step8:  Normalization
There is no special CPD class, all CPDs are implemented as Factors. This is intentional
Step9:  Marginalization
Sums out one or several variables. 
Again, note that the result is not necessarily normalized.
Step10:  Factor division
Step11:  Random sampling
Sampling from a Factor
Step12:  Sampling from DGM
Step13:  Fitting
DGMs
Step14:  Distributions in this DGM must be close to the real distributions (those in dgm).
Step15:  Factors
Step16:  Inference
Step17:  Structure Learning
Currently only Chow-Liu algorithm is implemented. It allows to find a tree structure with maximum likelihood.
Step18:  You can see that the Chow-Liu algorithm didn't find the true structure (since it is not a tree). However, the result is quite close to it.
Random DGM Generation
<ASSISTANT_TASK:>
Python Code:
dgm = gm.DGM.read('../networks/earthquake.bif')
dgm.draw() # you can move the cursor on a node to see it's CPD
Explanation: Directed Graphical Models
DGM is a class implementing Directed Graphical Models.
Loading DGMs
You can load bayesian networks in .bif format. For example, here we load a simple 'Earthquake' network.
Some networks can be downloaded here: http://www.bnlearn.com/bnrepository/
End of explanation
nx.dag_longest_path(dgm)
nx.average_neighbor_degree(dgm)
Explanation: DGM is a subclass of networkx.DiGraph, so you can use any networkx functions on it.
End of explanation
list(dgm.immoralities) # list of all immoralities in graph
list(dgm.v_structures) # list of all v_structures in graph
changed_dgm1 = gm.DGM()
changed_dgm1.add_edges_from(dgm.edges())
changed_dgm1.add_edge('JohnCalls', 'MaryCalls')
changed_dgm1.draw()
Explanation: Also, some DGM-specific queries about the graph are implemented:
End of explanation
changed_dgm1.is_I_equivalent(dgm)
changed_dgm1.is_I_map(dgm)
dgm.is_I_map(changed_dgm1)
Explanation: Checks for I-equivalence and I-maps
End of explanation
dgm.reachable('JohnCalls', observed=[])
dgm.reachable('JohnCalls', observed=['Alarm'])
Explanation: Find which nodes are reachable in the sense of d-separation
End of explanation
dgm.factor('JohnCalls') # get a factor assigned to variable JohnCalls
dgm.cpd('Alarm') # get a cpd of variable Alarm (i.e. normalized factor, see below)
Explanation: Factors
Currently, only TableFactor class is implemented. TableFactor represens an unnormalized discrete CPD.
End of explanation
alarm = dgm.factor('Alarm')
john_calls = dgm.factor('JohnCalls')
burglary = dgm.factor('Burglary')
print(alarm.arguments) # all variables in the network
print(alarm.scope) # variables in the factor
Explanation: Operations
All factor operations are either vectorized or written in C++, so they should be quite fast.
End of explanation
alarm * john_calls # note that the result is NOT a normalized CPD
Explanation: Factor Multiplication
End of explanation
(alarm * john_calls).normalize('JohnCalls') # now the result is a valid CPD: P(JohnCalls | Burglary, Earthquake, Alarm)
Explanation: Normalization
There is no special CPD class, all CPDs are implemented as Factors. This is intentional: often there is no need to normalize probability distributions on intermediate steps of an algorithm; doing that would hurt performance. 
.normalize([variables, ]) function normalized the factor so that the probabilities over variables sum to 1.
End of explanation
factor = burglary * alarm
factor
factor.marginalize('Burglary')
Explanation: Marginalization
Sums out one or several variables. 
Again, note that the result is not necessarily normalized.
End of explanation
alarm
factor = burglary * alarm
factor
factor / burglary
Explanation: Factor division
End of explanation
f = gm.TableFactor(arguments=['a', 'b'], scope=['a', 'b'])
f.table = np.array([[0.5, 0.5], [0.8, 0.2]])
f.normalize(*f.scope, copy=False)
f
# sampling
rv = f.rvs(size=1e+6)
# Check if the sampling is correct
f1 = gm.TableFactor(['a', 'b'], ['a', 'b']).fit(rv)
f1
# conditional sampling
f.rvs(observed={'a': [0, 1, 0]})
dgm.factor('Alarm').rvs(observed={
        'Earthquake': ['True', 'False', 'True'], 
        'Burglary': ['False', 'False', 'True']
    })
Explanation: Random sampling
Sampling from a Factor
End of explanation
dgm.draw()
dgm.rvs(size=10)
Explanation: Sampling from DGM
End of explanation
# generating data
data = dgm.rvs(size=1000000)
# copying graph
dgm1 = gm.DGM()
dgm1.add_nodes_from(dgm.nodes(), cpd=gm.TableFactor) # here we also need to specify cpd type (TableFactor)
dgm1.add_edges_from(dgm.edges())
dgm1.draw()
# fitting
dgm1.fit(data)
dgm1.draw()
Explanation: Fitting
DGMs
End of explanation
dgm.draw()
Explanation: Distributions in this DGM must be close to the real distributions (those in dgm).
End of explanation
factor = dgm.factor('Alarm') * dgm.factor('Burglary') * dgm.factor('Earthquake')
# generating data
data = factor.rvs(size=1000000)
new_factor = gm.TableFactor(['Alarm', 'Burglary', 'Earthquake'], ['Alarm', 'Burglary', 'Earthquake'])
# fitting
new_factor.fit(data)
factor
Explanation: Factors
End of explanation
# creating inference object
inference = gm.SumProductInference(dgm)
# making inference queries
inference(query=['JohnCalls']) # probability of John calling
# probability of John calling in case there was a burglary
inference(query=['JohnCalls'], observed=dict(Burglary='True'))
# probability of John/Mary calling in case there was a burglary
inference(query=['JohnCalls', 'MaryCalls'], observed=dict(Burglary='True'))
# probability of Earthquake given that both John and Mary called
inference(query=['Earthquake'], observed=dict(JohnCalls='True', MaryCalls='True'))
Explanation: Inference
End of explanation
# generating data
data = dgm.rvs(size=1000000)
learned = gm.chow_liu(data)
learned.draw()
Explanation: Structure Learning
Currently only Chow-Liu algorithm is implemented. It allows to find a tree structure with maximum likelihood.
End of explanation
gen = gm.ErdosRenyiDGMGen(factor_gen=gm.DirichletTableFactorGen())
gen().draw()
gen = gm.TreeDGMGen(factor_gen=gm.DirichletTableFactorGen())
gen().draw()
Explanation: You can see that the Chow-Liu algorithm didn't find the true structure (since it is not a tree). However, the result is quite close to it.
Random DGM Generation
End of explanation
<END_TASK>