<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Variant effect prediction
The variant effect prediction parts integrated in concise are designed to extract importance scores for a single nucleotide variant in a given sequence. Predictions are made for each output individually for a multi-task model. In this short tutorial we will be using a small model to explain the basic functionality and outputs.
At the moment there are three different effect scores to be chosen from. All of them require as in input
Step1:  As with any prediction that you want to make with a model it is necessary that the input sequences have to fit the input dimensions of your model, in this case the reference and alternative sequences in their forward and reverse-complement state have to have the shape [?, 1000, 4].
We will be storing the dataset in a dictionary for convenience
Step2:  All prediction functions have the same general set of required input values. Before going into more detail of the individual prediction functions We will look into how to run them. The following input arguments are availble for all functions
Step3:  The output of all functions is a dictionary, please refer to the individual chapters further on for an explanation of the individual values. Every dictionary contains pandas dataframes as values. Every column of the dataframe is named according to the values given in the out_annotation_all_outputs labels and contains the respective predicted scores.
Convenience function
For convenience there is also a function available which enables the execution of all functions in one call.
Additional arguments of the effect_from_model function are
Step4:  Again the returned value is a dictionary containing the results of the individual calculations, the keys are the names of the executed functions
<ASSISTANT_TASK:>
Python Code:
from effect_demo_setup import *
from concise.models import single_layer_pos_effect as concise_model
import numpy as np
# Generate training data for the model, use a 1000bp sequence
param, X_feat, X_seq, y, id_vec = load_example_data(trim_seq_len = 1000)
# Generate the model
dc = concise_model(pooling_layer="sum",
                   init_motifs=["TGCGAT", "TATTTAT"],
                   n_splines=10,
                   n_covariates=0,
                   seq_length=X_seq.shape[1],
                   **param)
# Train the model
dc.fit([X_seq], y, epochs=1,
       validation_data=([X_seq], y))
    
# In order to select the right output of a potential multitask model we have to generate a list of output labels, which will be used alongside the model itself.
model_output_annotation = np.array(["output_1"])
Explanation: Variant effect prediction
The variant effect prediction parts integrated in concise are designed to extract importance scores for a single nucleotide variant in a given sequence. Predictions are made for each output individually for a multi-task model. In this short tutorial we will be using a small model to explain the basic functionality and outputs.
At the moment there are three different effect scores to be chosen from. All of them require as in input:
The input sequence with the variant with its reference genotype
The input sequence with the variant with its alternative genotype
Both aformentioned sequences in reverse-complement
Information on where (which basepair, 0-based) the mutation is placed in the forward sequences
The following variant scores are available:
In-silico mutagenesis (ISM):
Predict the outputs of the sequences containing the reference and alternative genotype of the variant and use the differential output as a effect score.
Gradient-based score
Dropout-based score
Calculating effect scores
Firstly we will need to have a trained model and a set of input sequences containing the variants we want to look at. For this tutorial we will be using a small model:
End of explanation
import h5py
dataset_path = "%s/data/sample_hqtl_res.hdf5"%concise_demo_data_path
dataset = {}
with h5py.File(dataset_path, "r") as ifh:
    ref = ifh["test_in_ref"].value
    alt = ifh["test_in_alt"].value
    dirs = ifh["test_out"]["seq_direction"].value
    
    # This datset is stored with forward and reverse-complement sequences in an interlaced manner
    assert(dirs[0] == b"fwd")
    dataset["ref"] = ref[::2,...]
    dataset["alt"] = alt[::2,...]
    dataset["ref_rc"] = ref[1::2,...]
    dataset["alt_rc"] = alt[1::2,...]
    dataset["y"] = ifh["test_out"]["type"].value[::2]
    
    # The sequence is centered around the mutatiom with the mutation occuring on position when looking at forward sequences
    dataset["mutation_position"] = np.array([500]*dataset["ref"].shape[0])
Explanation: As with any prediction that you want to make with a model it is necessary that the input sequences have to fit the input dimensions of your model, in this case the reference and alternative sequences in their forward and reverse-complement state have to have the shape [?, 1000, 4].
We will be storing the dataset in a dictionary for convenience:
End of explanation
from concise.effects.ism import ism
from concise.effects.gradient import gradient_pred
from concise.effects.dropout import dropout_pred
ism_result = ism(model = dc, 
                 ref = dataset["ref"], 
                 ref_rc = dataset["ref_rc"], 
                 alt = dataset["alt"], 
                 alt_rc = dataset["alt_rc"], 
                 mutation_positions = dataset["mutation_position"], 
                 out_annotation_all_outputs = model_output_annotation, diff_type = "diff")
gradient_result = gradient_pred(model = dc, 
                                ref = dataset["ref"], 
                                ref_rc = dataset["ref_rc"], 
                                alt = dataset["alt"], 
                                alt_rc = dataset["alt_rc"], 
                                mutation_positions = dataset["mutation_position"], 
                                out_annotation_all_outputs = model_output_annotation)
dropout_result = dropout_pred(model = dc, 
                              ref = dataset["ref"], 
                              ref_rc = dataset["ref_rc"], 
                              alt = dataset["alt"], alt_rc = dataset["alt_rc"], mutation_positions = dataset["mutation_position"], out_annotation_all_outputs = model_output_annotation)
gradient_result
Explanation: All prediction functions have the same general set of required input values. Before going into more detail of the individual prediction functions We will look into how to run them. The following input arguments are availble for all functions:
model: Keras model
ref: Input sequence with the reference genotype in the mutation position
ref_rc: Reverse complement of the 'ref' argument
alt: Input sequence with the alternative genotype in the mutation position
alt_rc: Reverse complement of the 'alt' argument
mutation_positions: Position on which the mutation was placed in the forward sequences
out_annotation_all_outputs: Output labels of the model.
out_annotation: Select for which of the outputs (in case of a multi-task model) the predictions should be calculated.
The out_annotation argument is not required. We will now run the available predictions individually.
End of explanation
from concise.effects.snp_effects import effect_from_model
# Define the parameters:
params = {"methods": [gradient_pred, dropout_pred, ism],
         "model": dc,
         "ref": dataset["ref"],
         "ref_rc": dataset["ref_rc"],
         "alt": dataset["alt"],
         "alt_rc": dataset["alt_rc"],
         "mutation_positions": dataset["mutation_position"],
         "extra_args": [None, {"dropout_iterations": 60},
         		{"rc_handling" : "maximum", "diff_type":"diff"}],
         "out_annotation_all_outputs": model_output_annotation,
         }
results = effect_from_model(**params)
Explanation: The output of all functions is a dictionary, please refer to the individual chapters further on for an explanation of the individual values. Every dictionary contains pandas dataframes as values. Every column of the dataframe is named according to the values given in the out_annotation_all_outputs labels and contains the respective predicted scores.
Convenience function
For convenience there is also a function available which enables the execution of all functions in one call.
Additional arguments of the effect_from_model function are:
methods: A list of prediction functions to be executed. Using the same function more often than once (even with different parameters) will overwrite the results of the previous calculation of that function.
extra_args: None or a list of the same length as 'methods'. The elements of the list are dictionaries with additional arguments that should be passed on to the respective functions in 'methods'. Arguments defined here will overwrite arguments that are passed to all methods.
**argv: Additional arguments to be passed on to all methods, e.g,: out_annotation.
End of explanation
print(results.keys())
Explanation: Again the returned value is a dictionary containing the results of the individual calculations, the keys are the names of the executed functions:
End of explanation
<END_TASK>