<jupyter_start><jupyter_text>Generating ModelsThis notebook demonstrates how to generate a model from a compact, string representation of the transcription factors for a gene. We refer to this as a *gene descriptor* or just *descriptor*.There are 3 variants of descriptors, depending on the number of transcription factors:- 0 TF: n. Example: "1". This is just the gene number. A "0" is used for INPUT.- 1 TF: nsm. Example: "1+2". This indicates that P2 activates gene 1. "s" is sign, either "+" or "-".- 2 TF: nsmisk. Example: "6+7O-1". This means that P7 activates Gene 6 or Gene 1 inhibits Gene 6. The "O" is a non-competitve OR integration term (or "i") and the "1" corresponds to "k" in the specification of this descriptor. Other integration terms are: "A" for AND and "P" for competitive OR integration.<jupyter_code>%matplotlib inline
import tellurium as te
import gene_network as gn
import modeling_game as mg
import model_fitting as mf

import copy
import lmfit
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd<jupyter_output><empty_output><jupyter_text>Creating Gene ReactionsThis section is intended to give insight into how the code works. It is not essential for using the codes.Network generation is done by constructing from each gene descriptor a reaction for the mRNA and a set of constants, some of which may need to be initialize.<jupyter_code># Create the reaction (including kinetics) for gene7
gene7 = gn.GeneReaction(7)
print(gene7)
result = gn.GeneReaction.do("7+4")
result.mrna_reaction
result.constants
result = gn.GeneReaction.do("7+4A-1")
result.mrna_reaction<jupyter_output><empty_output><jupyter_text>Creating and Using a Generated ModelBelow are examples of running simulations using a network expressed in terms of gene descriptors, and a comparison between the simulation outputs of this representation and the original gene network.<jupyter_code># Descriptors for the initial network
print(gn.INITIAL_NETWORK)
# Construct the network. Starts with INITIAL_NETWORK
network = gn.GeneNetwork()
# Generate the network
network.generate()
# network.generate() sets two properties in network. The first is the Antimony model.
print(network.model)
# The second is the lmfit.Parameters() needed for parameter optimization
network.parameters
# Here's the result of simulating the generated network
rr = te.loada(network.model)
data = rr.simulate(0, 1200, 120)
df_data = mg.makeDF(data)
_ = plt.plot(df_data.index, df_data)
_ = plt.legend(df_data.columns)
def simulateInitialModel():
    """
    Simulates the original model. Returns a dataframe of the result with just mRNA.
    """
    with open("initial-Antimony-model.txt", "r") as fd:
        initial_model = fd.readlines()
    model = "".join(initial_model)
    rr = te.loada(model)
    data = rr.simulate(0, 1200, 120)
    return mg.makeDF(data)
# Compare the above results with the file initial-Antimony-model.txt
df_initial = simulateInitialModel()
_ = plt.plot(df_initial.index, df_initial)
_ = plt.legend(df_data.columns)
df = df_initial - df_data
df.sum().sum()<jupyter_output><empty_output><jupyter_text>Using Generated Networks in Network Discovery<jupyter_code># Observed data
# mRNA data for wild type
df_rnaseq = mf.cleanColumns(pd.read_csv("wild.csv"))
df_rnaseq.head()
# Here's the result of simulating the generated networ
mg.plotSimulation(df_rnaseq, network.model)
# Default initial network
gn.INITIAL_NETWORK
# Evaluating the initial model
network = gn.GeneNetwork()  # By default, uses the correct initial network
# Generate the network
network.generate()
# parameters = mg.runExperiment(df_rnaseq, network.model, network.parameters)
parameters = mg.runExperiment(df_rnaseq, network.model, ["Vm1"])
def evaluate(network, parameters, descriptors, sim_time=300, **kwargs):
    """
    Generates plots to evaluate the network using the parameters and descriptors.
    Only estimates the parameters for the new descriptors.
    :param GeneNetwork network:
    :param lmfit.Parameters parameters: parameters to use in constructing initializations
    :param list-str descriptors:
    :param int sim_time: length of the simululation
    :param dict kwargs: optional arguments for runExperiment
    :return GeneNetwork, lmfit.parameters:
    """
    new_network = network.copy()
    new_network.addInitialization(parameters)
    new_network.update(descriptors)
    # Generate the network
    new_network.generate()
    if len(new_network.new_parameters) == 0:
        news = "Vm1"
    else:
        news = new_network.new_parameters
    new_parameters = mg.makeParameters(news)
    # Just optimize the parameters that are initialized to 0, new_parameters
    newer_parameters = mg.runExperiment(df_rnaseq, new_network.model, new_parameters, 
                                      sim_time=sim_time, **kwargs)
    return new_network, newer_parameters

def compareModels(network, parameters, descriptors, sim_time=300):
    """
    Compares changes to a single gene. The typical use is to evaluate different
    regulations of a single gene.
    :param GeneNetwork network:
    :param lmfit.Parameters parameters: parameters to use in constructing initializations
    :param list-str descriptors: Gene descriptors for which plots are done
    :param int sim_time: length of the simululation
    """
    is_first = True
    for descriptor in descriptors:
        if is_first:
            is_plot_observations = True
            is_first = False
        else:
            is_plot_observations = False
        evaluate(network, parameters, [descriptor], sim_time=sim_time,
                is_plot_model=False, is_plot_observations=is_plot_observations,
                title=str(descriptor))
network = gn.GeneNetwork()  # By default, uses the correct initial network
network.update(['4-3','5+7'])
# Generate the network
network.generate()
# Just optimize the parameters that are initialized to 0, new_parameters
parameters1 = mg.runExperiment(df_rnaseq, network.model, network.new_parameters, 
                                  sim_time=300)
# network2, parameters2 = evaluate(network, parameters1, ['7-4'], sim_time=300)
compareModels(network, parameters1, ['7-7A-2'], sim_time=300)

new_network = network.copy()
new_network.addInitialization(parameters1)
new_network.update(['7-7A+4'])
new_network.generate()
# Just optimize the parameters that are initialized to 0, new_parameters
parameters2 = mg.runExperiment(df_rnaseq, new_network.model, new_network.new_parameters, 
                                  sim_time=300)
print(new_network.model)<jupyter_output>// Created by libAntimony v3.9.4

model *pathway()



  // Compartments and Species:

  species INPUT, P1, mRNA1, P2, mRNA2, P3, mRNA3, P4, mRNA4, P5, mRNA5, P6;

  species mRNA6, P7, mRNA7, P8, mRNA8;
J1: => mRNA1; L1 - d_mRNA1*mRNA1 + Vm1 * ( K1_1*INPUT^H1 + K2_1*P4^H1 + K3_1*INPUT^H1*P4^H1 ) / ( 1 + K1_1*INPUT^H1 + K2_1*P4^H1 + K3_1*INPUT^H1*P4^H1 )
J2: => mRNA2; L2 - d_mRNA2*mRNA2 + Vm2 * ( K1_2*P4^H2 ) / ( 1 + K1_2*P4^H2 )
J3: => mRNA3; L3 - d_mRNA3*mRNA3 + Vm3 * ( K1_3*P6^H3 ) / ( 1 + K1_3*P6^H3 )
J4: => mRNA4; L4 - d_mRNA4*mRNA4 + Vm4 * ( 1 ) / ( 1 + K1_4*P3^H4 )
J5: => mRNA5; L5 - d_mRNA5*mRNA5 + Vm5 * ( K1_5*P7^H5 ) / ( 1 + K1_5*P7^H5 )
J6: => mRNA6; L6 - d_mRNA6*mRNA6 + Vm6 * ( K1_6*P7^H6 ) / ( 1 + K1_6*P7^H6 + K2_6*P1^H6 + K3_6*P7^H6*P1^H6 )
J7: => mRNA7; L7 - d_mRNA7*mRNA7 + Vm7 * ( K2_7*P4^H7 ) / ( 1 + K1_7*P7^H7 + K2_7*P4^H7 + K3_7*P7^H7*P4^H7 )
J8: => mRNA8; L8 - d_mRNA8*mRNA8 + Vm8 * ( 1 ) / ( 1 + K1_8*P1^H8 )

  

// Protein Reactions:

  F1:  => P1; a_protein1*mRNA1[...]
{'fraction_non_alphanumeric': 0.0991811938721606, 'fraction_numerical': 0.03420496566296883, 'mean_word_length': 4.127285037237644, 'pattern_counts': {'":': 0, '<': 16, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '79639', 'license': 'MIT', 'n_tokens_mistral': 2650, 'n_tokens_neox': 2469, 'n_words': 759, 'path': 'archived_lectures/Fall_2019/lecture_17/lecture_17_Generating_Models.ipynb', 'repo_name': 'BioModelTools/topics-course', 'size': '7572'}