Hw18

README.md

@@ -1,163 +1,186 @@
+# BioProg_SV
 
-# BIOPROG_SV
+**BioProg_SV** (Biological Program by Stepnova Valeriya) is a toolkit designed for performing basic DNA/RNA sequence analysis and filtering FASTQ data based on quality metrics. The toolset includes functions for sequence transformation, filtering, and analysis.
 
+---
 
-**BioProg_SV** (Biological Programm by Stepnova Valeriya) is a toolkit designed to calculate basic DNA/RNA analysis and validate DNA quality from FastQC data.
-
-Authors:
-**Stepanova Valeria - Bioinformatics institute 2024/2025**
-
+## Authors
 
-**BioProg_SV** documentation is available on the GitHub [repo](https://github.com/Stepanovalera/BioProgSV).<br/>
+- **Stepanova Valeriya**  
+  Bioinformatics Institute, 2024–2025
 
+Documentation is available in the GitHub repository: [Repository Link](https://github.com/your-repo-link).
 
-## Content
+---
 
+## Table of Contents
 
-* [Instructions](#instructions)
-* [Example input](#examples)
-* [Contact](#contact)
+1. [Instructions](#instructions)
+2. [Function Descriptions](#function-descriptions)
+3. [Usage Examples](#usage-examples)
+4. [Contact](#contact)
 
+---
 
 ## Instructions
-### Function Descriptions
-This program contains main functions:
 
-#### `filter_fastqc()`
-This function processes FastQC output data to filter and validate DNA sequences based on quality metrics. It applies multiple quality criteria to determine which sequences from a given dictionary (seqs) pass the filtering steps. The filtering criteria include:
+### Installing Dependencies
 
-**Description:**
-Filters sequences in a FASTQ file based on GC content, length, and quality threshold.
+If the project uses third-party libraries, install them using the `requirements.txt` file:
 
-**Parameters:**
+```bash
+pip install -r requirements.txt
+```
 
-- `input_fastq` (str): Path to the input FASTQ file.
-
-- `output_fastq` (str): Name of output FASTQ file.
-
-- `gc_bounds` (tuple or int): GC content bounds for filtering.
+If the file is empty, it means the project only uses Python's standard libraries.
 
-- `length_bounds` (tuple or int): Length bounds for filtering.
+### Running the Program
 
-- `quality_threshold`(float): Minimum average quality score for filtering.
+The program supports two main modes of operation:
 
-  *GC Content Bounds*: Only sequences within the specified percentage range are retained. The GC content influences sequence stability and should fall within acceptable bounds for specific analyses.  Default paramenter: `gc_bounds = (0, 100)` Input can contain a single number, which will be interpreted as the upper bound. 
+1. **DNA/RNA Sequence Analysis**:
+   ```bash
+   python main.py dna_rna_tools --seqs "ATG" --action transcribe
+   ```
 
-  *Length Bounds*: The function checks if sequence lengths are within the specified minimum and maximum values, ensuring that only sequences of appropriate length are included. Default paramenter: `length_bounds = (0, 2**32)` Input can contain a single number, which will be interpreted as the upper bound.
+2. **FASTQ File Filtering**:
+   ```bash
+   python main.py filter_fastq input.fastq output.fastq --gc_bounds 50 60 --length_bounds 50 100 --quality_threshold 30
+   ```
 
-  *Quality Threshold*: Sequences must exceed a minimum quality score, ensuring that only high-quality data are used for downstream analyses. Default paramenter: `quality_threshold = 0`
+---
 
-**Returns:**
+## Function Descriptions
 
-output FASTA file
+### 1. `filter_fastq`
 
+This function filters sequences from a FASTQ file based on specified quality criteria.
 
+#### Parameters:
+- `input_fastq` (str): Path to the input FASTQ file.
+- `output_fastq` (str): Path to the output FASTQ file.
+- `gc_bounds` (tuple or int): GC content bounds (default: `(0, 100)`). If a single number is provided, it is treated as the upper bound.
+- `length_bounds` (tuple or int): Sequence length bounds (default: `(0, 2**32)`). If a single number is provided, it is treated as the upper bound.
+- `quality_threshold` (float): Minimum average quality score required for sequences to be included (default: `0.0`).
 
-  in-build function `fast_qc` - The function calculates:
+#### Returns:
+- A filtered FASTQ file containing sequences that meet the specified criteria.
 
+---
 
-    - GC Content: Percentage of guanine (G) and cytosine (C) in the DNA sequence.
+### 2. `fast_qc`
 
+This function calculates quality metrics for DNA/RNA sequences.
 
-    - Sequence Length: The number of bases in the DNA sequence.
+#### Metrics Calculated:
+- **GC Content**: Percentage of guanine (G) and cytosine (C) in the sequence.
+- **Sequence Length**: Number of bases in the sequence.
+- **Average Quality Score**: Mean quality score derived from Phred quality scores.
 
+#### Parameters:
+- `seqs` (dict): A dictionary where keys are sequence names, and values are tuples containing the sequence (str) and its quality string (str).
 
-    - Average Quality Score: An average quality derived from Phred quality scores
+#### Returns:
+- A dictionary with sequence names as keys and tuples as values. Each tuple contains:
+  - GC content (float)
+  - Sequence length (int)
+  - Average quality score (float)
 
+---
 
-#### `convert_multiline_fasta_to_oneline(input_fasta, output_fasta)`
+### 3. `convert_multiline_fasta_to_oneline`
 
-**Description:**  
-This function processes a FASTA file, converting sequences that are split over multiple lines into a single line format. Each sequence starts with a header line beginning with `>`. The function ensures that every sequence is represented as a single line in the output file.
+Converts a multi-line FASTA file into a single-line format.
 
-**Parameters:**
-- `input_fasta` (str): The path to the input FASTA file, where sequences may be spread across multiple lines.
-- `output_fasta` (str): The reformatted FASTA file name.
+#### Parameters:
+- `input_fasta` (str): Path to the input FASTA file.
+- `output_fasta` (str): Path to the output FASTA file.
 
-**Returns:**
-- (str):  output FASTA file where sequences are converted to a one-line format.
+#### Returns:
+- A FASTA file where each sequence is represented as a single line.
 
 ---
 
-#### `parse_blast_output(input_blast, output_blast)`
+### 4. `parse_blast_output`
 
-**Description:**  
-This function parses the output of a BLAST search, extracting relevant sequence descriptions. It starts recording after encountering the 'Description' line and continues until it reaches the first empty line. The function then saves a truncated version of each description line, containing only the portion before the ellipsis (`...`).
+Parses BLAST output files to extract relevant sequence descriptions.
 
-**Parameters:**
-- `input_blast` (str): The path to the input BLAST output file, which contains the results to be parsed.
-- `output_blast` (str): The file name where the parsed BLAST output will be saved.
+#### Parameters:
+- `input_blast` (str): Path to the input BLAST output file.
+- `output_blast` (str): Path to the output file for parsed results.
 
-**Returns:**
-- (str): The path to the output file containing the parsed BLAST results.
+#### Returns:
+- A file containing truncated descriptions of sequences.
 
 ---
 
-#### `select_genes_from_gbk_to_fasta(input_gbk, output_fasta, genes, n_before=1, n_after=1)`
+### 5. `select_genes_from_gbk_to_fasta`
 
-**Description:**  
-This function extracts protein-coding sequences from a GenBank (GBK) file and writes neighboring genes' sequences to a FASTA file, excluding the target genes themselves. It identifies genes of interest and writes their adjacent genes' sequences, defined by `n_before` and `n_after`, to the output file.
+Extracts neighboring gene sequences from a GenBank (GBK) file.
 
-**Parameters:**
-- `input_gbk` (str): The path to the input GenBank file, containing the sequence data to be processed.
-- `output_fasta` (str):  FASTA file name where selected genes' sequences will be saved.
-- `genes` (any): A list of target gene names adjacent to which the sequences will be extracted.
-- `n_before` (int, optional): The number of genes before each target gene to include in the output. Defaults to 1.
-- `n_after` (int, optional): The number of genes after each target gene to include in the output. Defaults to 1.
+#### Parameters:
+- `input_gbk` (str): Path to the input GenBank file.
+- `output_fasta` (str): Path to the output FASTA file.
+- `genes` (list): List of target gene names.
+- `n_before` (int, optional): Number of genes before the target gene to include (default: `1`).
+- `n_after` (int, optional): Number of genes after the target gene to include (default: `1`).
 
-**Returns:**
-- (str): the output FASTA file containing sequences of neighboring genes.
+#### Returns:
+- A FASTA file containing sequences of neighboring genes.
 
+---
 
+### 6. `run_dna_rna_tools`
 
-#### `run_dna_rna_tools() `
+Performs various DNA/RNA operations, including transcription, reverse complementation, and more.
 
-This function perform multiple DNA/RNA operations, providing a streamlined interface to process sequences efficiently.The function is case-insensitive.
+#### Supported Actions:
+- `transcribe`: Converts DNA to RNA.
+- `reverse`: Reverses a sequence.
+- `complement`: Generates the complement of a DNA sequence.
+- `reverse_complement`: Generates the reverse complement of a DNA sequence.
 
-  `is_none()` - Checks if the given sequence is empty or None. This utility function ensures that subsequent operations are performed only on valid sequences, preventing errors in data processing.
-
-  `is_DNA()` - Validates whether the provided sequence is a DNA sequence, checking for the presence of valid DNA nucleotides (A, T, C, G).
-
-  `is_RNA()` - Determines if the sequence is an RNA sequence by verifying the presence of valid RNA nucleotides (A, U, C, G). 
+#### Parameters:
+- `seqs` (list): List of DNA/RNA sequences.
+- `action` (str): Action to perform (e.g., `"transcribe"`, `"reverse"`).
 
-  `reverse()` - Reverses a given nucleotide sequence.
-  
+#### Returns:
+- A single sequence (if one input) or a list of processed sequences.
 
-  `reverse_complement()` - Generates the reverse complement of a DNA sequence.
-
+---
 
-  `transcribe()` - Converts a DNA sequence into an RNA sequence by replacing thymine with uracil.
-
-
-  `complement()` - Produces the complement of a DNA strand by replacing each nucleotide with its complementary base.
-
+## Usage Examples
 
-Each sequence is processed, and results are collected into a list. If there is only one sequence to process, the function directly returns the single result; otherwise, it returns a `list()`. 
+### Example 1: DNA/RNA Tools
 
-## Examples
+```python
+# Transcribe DNA to RNA
+run_dna_rna_tools(["ATG"], "transcribe")  # Output: ["AUG"]
 
-For `run_dna_rna_tools()`:
+# Reverse a sequence
+run_dna_rna_tools(["AGT"], "reverse")  # Output: ["TGA"]
 
+# Generate reverse complement
+run_dna_rna_tools(["GTGT"], "reverse_complement")  # Output: ["ACAC"]
+```
 
-~~~
-    run_dna_rna_tools("ATG", "transcribe") == "AUG"
-    run_dna_rna_tools("AGt", "ATTCC" "reverse") == "AGu"
-    run_dna_rna_tools("GTGT", "TGU" "complement") == "AGu"
-~~~
+### Example 2: FASTQ Filtering
 
-For `bioprogsv()`:
+```python
+# Filter sequences based on length
+filter_fastq("input.fastq", "output.fastq", length_bounds=(50, 100))
 
+# Filter sequences based on GC content and quality
+filter_fastq("input.fastq", "output.fastq", gc_bounds=(50, 60), quality_threshold=30)
+```
 
-~~~
-EXAMPLE_FASTQ = {
-    # 'name' : ('sequence', 'quality')
-    '@SRX079801': ('ACAGCAACATAAACATGATGGGATGGCGTAAGCCCCCGAGATATCAGTTTACCCAGGATAAGAGATTAAATTATGAGCAACATTATTAA', 'FGGGFGGGFGGGFGDFGCEBB@CCDFDDFFFFBFFGFGEFDFFFF;D@DD>C@DDGGGDFGDGG?GFGFEGFGGEF@FDGGGFGFBGGD'),
-    '@SRX079802': ('ATTAGCGAGGAGGAGTGCTGAGAAGATGTCGCCTACGCCGTTGAAATTCCCTTCAATCAGGGGGTACTGGAGGATACGAGTTTGTGTG', 'BFFFFFFFB@B@A<@D>BDDACDDDEBEDEFFFBFFFEFFDFFF=CC@DDFD8FFFFFFF8/+.2,@7<<:?B/:<><-><@.A*C>D'),
-    '@SRX079803': ('GAACGACAGCAGCTCCTGCATAACCGCGTCCTTCTTCTTTAGCGTTGTGCAAAGCATGTTTTGTATTACGGGCATCTCGAGCGAATC', 'DFFFEGDGGGGFGGEDCCDCEFFFFCCCCCB>CEBFGFBGGG?DE=:6@=>A<A>D?D8DCEE:>EEABE5D@5:DDCA;EEE-DCD')}
+---
 
-  filter_fastq(EXAMPLE_FASTQ, length_bounds=1)
-  filter_fastq(EXAMPLE_FASTQ, gc_bounds = (55, 90), quality_threshold=10)
-~~~
 ## Contact
 
-Also, you can send your feedback to [ukrainskaya49@gmail.com](mailto:ukrainskaya49@gmail.com).
+For questions, feedback, or suggestions, please contact:
+
+- Email: ukrainskaya49@gmail.com
+
+---
+

addscript/fastq_qc_function.py → addscript/qc_tool/fastq_qc_function.py

@@ -25,10 +25,16 @@ def fast_qc(seqs):
       indicating the reliability of the sequence data.
     """
     gc_len_q = {}
-    quality_scores = []
     for sequence_name, (sequence, quality) in seqs.items():
-        gc_count = (sequence.count('G') + sequence.count('C'))/len(sequence) * 100
-        quality_scores = [ord(char) - 33 for char in quality]
-        average_quality = sum(quality_scores) / len(quality_scores) if quality_scores else 0
-        gc_len_q[sequence_name] = (gc_count, len(sequence), average_quality)
+        if not sequence:  # Проверка на пустую последовательность
+            gc_count = 0.0
+            length = 0
+            average_quality = 0.0
+        else:
+            gc_count = (sequence.count('G') + sequence.count('C')) / len(sequence) * 100
+            quality_scores = [ord(char) - 33 for char in quality]
+            average_quality = sum(quality_scores) / len(quality_scores) if quality_scores else 0
+            length = len(sequence)
+
+        gc_len_q[sequence_name] = (gc_count, length, average_quality)
     return gc_len_q

addscript/qc_tool/fastq_read_write_file.py

@@ -0,0 +1,28 @@
+import os
+
+def read_fastq(input_fastq):
+    base_directory = os.path.dirname(input_fastq)
+    filtered_directory = os.path.join(base_directory, 'filtered')
+    os.makedirs(filtered_directory, exist_ok=True)
+
+    input_fastq_data = {}
+    with open(input_fastq) as file:
+        while True:
+            header = file.readline().strip()
+            if not header:
+                break
+            sequence = file.readline().strip()
+            _ = file.readline().strip()  # Пропускаем строку '+'
+            quality = file.readline().strip()
+            input_fastq_data[header] = (sequence, quality)
+
+    return input_fastq_data, filtered_directory
+
+
+def write_fastq(output_fastq_data, output_fastq):
+    with open(output_fastq, 'w') as file:
+        for sequence_id, (sequence, quality) in output_fastq_data.items():
+            file.write(f"{sequence_id}\n")
+            file.write(f"{sequence}\n")
+            file.write(f"+{sequence_id[1:]}\n")
+            file.write(f"{quality}\n")