File input for BSAT

Bio_Files_Processor.py

@@ -0,0 +1,126 @@
+import os
+from typing import List
+
+def convert_multiline_fasta_to_oneline(input_fasta: str, output_fasta: str = None) -> str:
+    """
+    Function conver multiline fasta file into fasta with name line marked by '>' and the other line with sequence.
+    Results save into new folder "Converted_data". If the folder already exists, new data is written to it.
+
+    :param input_fasta: Path to fasta-file with your seqs with. 
+     !! It is necessary to indicate the name along with extensions (.fasta) !!
+    :type input_fasta: str
+    :param output_fasta: Name of new fasta-file with your seqs in one line, default value = None
+    It is necessary to indicate the name along with extensions (.fasta)
+    :type output_fasta: str
+    :rtype: str
+    :return: Script completion message 
+    """
+    if input_fasta.find('.fasta') == 0:
+        raise ValueError(f'Wrong file format in input!')
+    if os.path.exists(os.path.join('.', 'Converted_data')) == False:
+        os.mkdir(os.path.join('.', 'Converted_data'))
+    if output_fasta == None:
+        out_name = 'one_line_' + input_fasta.strip("/")[-1]
+    else:
+        out_name = output_fasta
+    gene = 0
+    seq = 0
+    gene_name = ''
+    prot_seq = ''
+    gene_and_seq = dict()
+    with open (input_fasta) as seq_fasta:
+        for line in seq_fasta:
+            if gene == 1:
+                if line.startswith('>') == False:
+                    gene_and_seq[gene_name] += line.strip('\n')
+                else:
+                    gene_name = line.strip('\n')
+                    gene_and_seq[gene_name] = ''
+            if gene == 0 and line.startswith('>'):
+                gene_name = line.strip('\n')
+                gene_and_seq[gene_name] = ''
+                gene = 1
+    with open (os.path.join('.', 'Converted_data', out_name), mode = 'w') as new_seq_fasta:
+        for key, item in gene_and_seq.items():
+            new_seq_fasta.write(key + '\n')
+            new_seq_fasta.write(item + '\n')
+    return 'All sequences processed!'
+
+
+def select_genes_from_gbk_to_fasta(input_gbk: str, genes: List[str], n_before: int = 1, n_after: int = 1, 
+                                   output_fasta: str = None) -> str:
+    '''
+    Function help to search neighbours of GOI (gene of interest). Function writes neighbours of GOI in new FASTA-file 
+    as: name of gene, protein sequence. Results save into new folder "Analyzed_data". If the folder already exists, 
+    new data is written to it.
+
+    :input_gbk: Path to gbk-file with your seqs with.
+    !! It is necessary to indicate the name along with extensions (.gbk) !!
+    :type input_gbk: str
+    :param genes: Gene of interest names
+    :type genes: List[str]
+    :param n_before: number of genes before GOI (>0), default value = 1
+    :type n_before: int
+    :param n_after: number of genes after GOI (>0), default value = 1
+    :type n_after: int
+    :output_fasta: Name of FASTA-file with neighbours of GOI (names and seqs), default value = None
+    :type output_fasta: str
+    :rtype: str
+    :return: Script completion message
+
+    '''
+    if input_gbk.find('.gbk') == 0:
+        raise ValueError(f'Wrong file format in input!')
+    if os.path.exists(os.path.join('.', 'Analyzed_data')) == False:
+        os.mkdir(os.path.join('.', 'Analyzed_data'))
+    genes_for_search = genes
+    genes_gbk = []
+    genes_for_search_in_gbk = []
+    neighbour_genes = dict()
+    if output_fasta == None:
+        output_fasta = 'output_for_gbk.fasta'
+    with open (input_gbk) as gbk:
+        for line in gbk:
+            if '/gene' in line:
+                genes_gbk += [line.strip().split('=')[1]]
+    for el in genes_for_search:
+        genes_for_search_in_gbk += [gn for gn in genes_gbk if el in gn]
+    for gene in genes_for_search_in_gbk:
+        gene_index = genes_gbk.index(gene)
+        if gene_index >= 0 and gene_index < (len(genes_gbk) - 1):
+            for i in range(1, n_before + 1):
+                neighbour_genes[(genes_gbk[gene_index - i])] = 0
+            for i in range(1, n_after + 1):
+                neighbour_genes[(genes_gbk[gene_index + i])] = 0
+        else:
+            for i in range(1, n_before + 1):
+                neighbour_genes[(genes_gbk[gene_index - i])] = 0
+            for i in range(1, n_after):
+                neighbour_genes[(genes_gbk[0 + i])] = 0
+    with open (input_gbk) as gbk:
+        gene_name_read = 0
+        protein_read = 0
+        gene_name = ''
+        protein = ''
+        for line in gbk:
+            if '/gene' in line:
+                gene_name = line.strip().split('=')[1]
+                if gene_name in neighbour_genes:
+                    gene_name_read = 1
+            if protein_read == 1:
+                protein += line.strip('\n').strip(' ')
+                if '"' in line:
+                    protein_read = 0
+                    gene_name_read = 0
+                    neighbour_genes[gene_name] = protein
+                    protein = ''
+            if gene_name_read == 1 and '/translation' in line:
+                protein_read = 1
+                protein += line.strip().split('=')[1]
+    with open (os.path.join('.', 'Analyzed_data', output_fasta), mode = 'w') as fasta:
+        for name, seq in neighbour_genes.items():
+            name = '>'+name.replace('\"','')
+            fasta.write(name + '\n')
+            fasta.write(seq.replace('\"','') + '\n')
+    return 'All sequences processed!'
+

Bio_Seq_Analysis_Tool.py

@@ -0,0 +1,163 @@
+from typing import Union, Tuple, List, Set, Dict
+from modules_for_BSAT import fastq_analysis as fq
+from modules_for_BSAT  import dna_rna_analysis as na
+from modules_for_BSAT import protein_analysis as pa
+
+
+def dna_rna_analysis(*args: str, operation: str) -> Union[List[float], List[str]]:
+    """
+    This function performs a number of operations on DNA or RNA.
+
+    Operations supported by this functions:
+    -transcribe - return transcribed sequence
+    -reverse - return reverse sequence
+    -complement - return complement sequence
+    -reverse_complement - return reverse complement sequence
+    -gc_calculate - return sequence GC-content in percent
+
+    :param args: nucleic acid sequence
+    :type param seqs: str
+
+    :param operation: type of operation required
+    :type param operation: str
+
+    :return: Analysis of nucleic acid sequence
+    :rtype: List[str]
+
+    :raises ValueError: if sequence not RNA or DNA, also if the operation value out of OPERATION_DICT
+    """
+    OPERATION_DICT = {'transcribe': na.transcribe, 
+                'reverse': na.reverse,
+                'reverse_complement': na.reverse_complement,
+                'complement': na.complement,
+                'gc_calculate': na.gc_calculate}
+    analysis = []
+    for seq in args:
+        na.is_na(seq)
+        if operation in OPERATION_DICT.keys():
+            analysis.append(OPERATION_DICT[operation](seq))
+        else:
+            raise ValueError(f'Wrong operation!')
+    return analysis
+
+
+def analyse_fastq(input_path: str, 
+                  gc_bounds: Union[int, float, Tuple [int], Tuple [float]] = (0, 100), 
+                  length_bounds: Union[int, Tuple [int]] = (0, 2**32),
+                  quality_threshold: float = 0.0, filtered_file_name: Union[None, str] = None) -> Dict[str,str]:
+    """
+    This function help analyze a set of reads obtained from next-generation sequencing. 
+
+    The function allow to filter the desired reads according to three parameters:
+    GC-content, length and reading quality.
+
+    :param seqs: 
+    Path to the file with FASTQ-sequences in the format. 
+    :type seqs: str
+
+    :param gc_bounds: 
+    Boundary parameters for filtering sequences by GC-content. Save only reads with a GC-content between boundaries 
+    or lower than one boundary. Lower boundary cannot be less than 0 and upper boundary cannot be greater than 100. 
+    gc_bounds default value is (0,100).
+    :type param gc_bounds: Union[int, float, Tuple [int], Tuple [float]]
+
+    :param length_bounds: 
+    Boundary parameters for filtering sequences by length. Works the same as gc_bounds. Lower boundary cannot be less 
+    than 0 and upper boundary cannot be greater than 2^32. length_bounds default value is (0,2^32)
+    :type param length_bounds: Union[int, Tuple [int]
+
+    :param quality_threshold: 
+    Threshold for quality of each nucleotide in read. Quality incodes by ASCII codes. The threshold cannot be more 
+    than 40. quality_threshold default value is 0 
+    :type param quality_threshold: float
+
+    :return: 
+    New dictionaries with fastq sequence.The first one consisting of filtered fastq sequences and the other one with 
+    sequences that did not pass filters.
+    :rtype: Dict[str]
+
+    :raises ValueError: if sequence not RNA or DNA, also if the argument values are outside the allowed ones
+    """
+    seqs, path_to_file = fq.read_fastq(input_path)
+    file_name = path_to_file.split("/")[-1]
+    if type(gc_bounds) == float or type(gc_bounds) == int:
+        gc_bounds = (0,gc_bounds)
+        if gc_bounds[0] < 0 or gc_bounds[1] > 100:
+            raise ValueError(f'Wrong boundaries!')
+    if type(length_bounds) == int:
+        length_bounds = (0,length_bounds)
+        if length_bounds[0] < 0 or length_bounds[1] > 2**32:
+            raise ValueError(f'Wrong boundaries!')        
+    if quality_threshold > 40:
+        raise ValueError(f'Wrong quality threshold!')
+    analysed_seq = {}
+    error_seq = {}
+    for seq in seqs.items():
+        if fq.analyse_gc(seq[1][0]) >= gc_bounds[0] and fq.analyse_gc(seq[1][0]) <= gc_bounds[1]:
+            if fq.analyse_length(seq[1][0]) >= length_bounds[0] and fq.analyse_length(seq[1][0]) <= length_bounds[1]:
+                if fq.analyse_quality(seq[1][1]) > quality_threshold:
+                    analysed_seq[seq[0]] = (seq[1])
+                else:
+                    error_seq[seq[0]] = (seq[1])
+            else:
+                error_seq[seq[0]] = (seq[1])
+        else:
+            error_seq[seq[0]] = (seq[1])
+    if filtered_file_name == None:
+        new_file_name = file_name 
+    else:
+        new_file_name = filtered_file_name
+    fq.write_fastq(analysed_seq, path_to_file, new_file_name)
+    return error_seq
+
+
+def run_protein_analysis(*args: str) -> Union[List[str], List[float]]:
+    """
+    Launch desired operation with proteins sequences. Pass comma-separated sequences,
+    additional argument (if certain function requires it) and specify function name you want to apply to all
+    sequences.
+    Pass arguments strictly in this order, otherwise it won't be parsed.
+    If the input is a sequence of amino acids written in a three-letter code, then the amino acids must be separated 
+    by a space. If the input is a sequence of amino acids written in a single-letter code, then the amino acids may 
+    not be separated by a space.
+
+    :param args:
+    - seq (str): amino acids sequences for analysis in 1-letter or 3-letter code (as many as you wish)
+    - additional arg (str): necessary parameter for certain functions (for example, specify target protein site)
+    - operation name (str): specify procedure you want to apply
+    :type param args: str
+
+    :return: the result of procedure in list or str format
+    :rtype: Union[List[str], List[float]
+    :raises ValueError: if sequence not protein or sequence in wrong format    
+    """
+    # first value is function name, second is real function, third is number of function arguments
+    function_names = {'change_residues_encoding': [pa.change_residues_encoding, 2],
+                      'is_protein': [pa.is_protein, 1],
+                      'get_seq_characteristic': [pa.get_seq_characteristic, 1],
+                      'find_residue': [pa.find_residue, 2],
+                      'find_site': [pa.find_site, 2],
+                      'calculate_protein_mass': [pa.calculate_protein_mass, 1],
+                      'calculate_average_hydrophobicity': [pa.calculate_average_hydrophobicity, 1],
+                      'get_mrna': [pa.get_mrna, 1],
+                      'calculate_isoelectric_point': [pa.calculate_isoelectric_point, 1],
+                      'analyze_secondary_structure': [pa.analyze_secondary_structure, 1]}
+
+    procedure = args[-1]
+
+    processed_result = []
+
+    seqs = [pa.change_residues_encoding(seq) for seq in args[:-1 * (function_names[procedure][1])]]
+    for idx, seq in enumerate(seqs):
+        if not pa.is_protein(seq):
+            raise ValueError(f'Sequence {seq} is not protein!')
+            continue
+        if function_names[procedure][1] == 1:
+            processed_result.append(function_names[procedure][0](seq))
+        elif function_names[procedure][1] == 2:
+            add_arg = args[-2]
+            processed_result.append(function_names[procedure][0](seq, add_arg))
+    if len(processed_result) == 1:
+        return processed_result[0]
+    return processed_result
+