Lesson · 40 min · Free
Python Interlude: Scripting for Bioinformatics
Python Interlude: Scripting for Bioinformatics body { font-family: sans-serif; line-height: 1.6; margin: 20px; } h1, h2 { color: #2c3e50; } pre { background-color: #ecf0f1; padding: 15px; border-radius: 5px; overflow-x:
Python Interlude: Scripting for Bioinformatics
Welcome to this interlude module on Python scripting, specifically tailored for bioinformatics applications. While our course emphasizes command-line tools and foundational computational concepts, Python serves as an indispensable glue language for automating tasks, processing data, and building more complex pipelines in biomedicine. For pharmacy and biotech students, understanding how to leverage Python for data manipulation, file parsing, and basic statistical analysis is a critical skill that bridges the gap between raw data and actionable insights. Python's readability, extensive libraries (like Biopython, NumPy, Pandas), and vibrant community make it an ideal choice for bioinformatics. You'll find it used for everything from parsing large genomic files (e.g., FASTA, FASTQ, VCF) to orchestrating complex workflows that involve multiple external programs. This module will not make you a Python expert, but it will equip you with the fundamental concepts and practical examples necessary to start writing your own scripts to streamline your computational tasks.
Basic File I/O and String Manipulation
A common task in bioinformatics is reading data from files, processing each line or record, and then writing results to a new file. Python makes this straightforward with its built-in file handling capabilities. Let's look at an example where we read a simulated FASTA file, extract sequence IDs, and calculate the length of each sequence. # example_fasta.py # Create a dummy FASTA file for demonstration with open("sequences.fasta", "w") as f: f.write(">Seq1_Human_Chromosome1\n") f.write("ATGCGTACGTAGCTAGCTAGCTAGCTACGTAGCTACGTAGCTACGTAGCTA\n") f.write(">Seq2_Mouse_GeneX\n") f.write("GCTAGCTAGCTAGCTACGTAGCTACGTAGCTACGTAGCTAGCTAGCTACGTAGCTACGTAGCTACGTAGCTA\n") f.write(">Seq3_Ecoli_Plasmid\n") f.write("CGTAGCTACGTAGCTACGTAGCTACGTAGCTACGTAGCTACGTAGCTACGTAGCTACGTAGCTACGTAGCTA\n") print("Processing sequences.fasta...") sequence_data = {} current_id = None with open("sequences.fasta", "r") as infile: for line in infile: line = line.strip() # Remove leading/trailing whitespace, including newline characters if line.startswith(">"): current_id = line[1:] # Remove the '>' character sequence_data[current_id] = "" elif current_id: sequence_data[current_id] += line print("\nSequence Information:") for seq_id, sequence in sequence_data.items(): print(f"ID: {seq_id}, Length: {len(sequence)}") # Write sequence IDs and lengths to a new CSV file with open("sequence_summary.csv", "w") as outfile: outfile.write("Sequence_ID,Length\n") for seq_id, sequence in sequence_data.items(): outfile.write(f"{seq_id},{len(sequence)}\n") print("\nSequence summary written to sequence_summary.csv") In this script, we first simulate creating a FASTA file. Then, we open it for reading. The for line in infile: loop reads the file line by line, which is memory-efficient for large files. We use line.strip() to clean up whitespace and startswith('>') to identify header lines. String slicing ( line[1:] ) is used to remove the '>' symbol. Finally, we iterate through our collected data and print it, then write it to a simple CSV file. Another common pattern is to process a list of files. Imagine you have multiple output files from a command-line tool, and you want to extract a specific piece of information from each. Python's os module is incredibly useful for interacting with the operating system, including listing directories and constructing file paths. # process_results.py import os # Create dummy result files for demonstration if not os.path.exists("results"): os.makedirs("results") with open("results/sample_A_output.txt", "w") as f: f.write("Processing complete.\n") f.write("Detected variants: 15\n") f.write("Average depth: 120x\n") with open("results/sample_B_output.txt", "w") as f: f.write("Processing complete.\n") f.write("Detected variants: 22\n") f.write("Average depth: 95x\n") print("Extracting variant counts from result files...") variant_counts = {} results_dir = "results" for filename in os.listdir(results_dir): if filename.endswith("_output.txt"): filepath = os.path.join(results_dir, filename) with open(filepath, "r") as f: for line in f: if "Detected variants:" in line: # Extract the number using string splitting count_str = line.split(":")[1].strip() sample_name = filename.replace("_output.txt", "") variant_counts[sample_name] = int(count_str) break # Assuming only one "Detected variants" line per file print("\nVariant Counts per Sample:") for sample, count in variant_counts.items(): print(f"{sample}: {count} variants") # Clean up dummy files and directory os.remove("results/sample_A_output.txt") os.remove("results/sample_B_output.txt") os.rmdir("results") print("\nCleaned up dummy result files and directory.") This second example demonstrates how to iterate through files in a directory using os.listdir() and os.path.join() to construct full file paths. It then opens each relevant file, searches for a specific string ("Detected variants:"), extracts the numerical value, and stores it in a dictionary. This pattern is extremely powerful for summarizing results from numerous experiments or analyses.
Key Takeaways
Python is a crucial scripting language for automating bioinformatics tasks and data manipulation. File I/O ( open() , read() , write() , with open(...) as f: ) is fundamental for handling biological data files. String methods like strip() , startswith() , split() , and string slicing are essential for parsing text-based data. The os module provides functions for interacting with the file system (e.g., os.listdir() , os.path.join() , os.makedirs() , os.remove() , os.rmdir() ). Dictionaries are excellent data structures for storing key-value pairs, such as sequence IDs and their data, or sample names and their statistics.
Practice Exercise: Parse a Simple CSV File
Imagine you have a CSV file named gene_expression.csv with the following content: GeneID,SampleA,SampleB,SampleC GAPDH,1200,1150,1300 ACTB,800,820,790 TP53,50,55,48 MYC,200,210,195 Your task is to write a Python script that: Creates this gene_expression.csv file programmatically (like in the examples). Reads the file line by line. Parses each line to extract the GeneID and its expression values across samples. Calculates the average expression for each gene across all samples. Prints the GeneID and its calculated average expression in a user-friendly format (e.g., "Gene GAPDH: Average Expression = 1216.67"). Cleans up the created gene_expression.csv file. This exercise will reinforce file I/O, string splitting, and basic data type conversion (strings to numbers).
Watch the full lesson — free
This topic is part of Computational Biomedicine: From Command Line to Single-Cell, a complete AI-narrated video course. Press play once and watch the entire lecture like a movie.
Start the course free →