|
| 1 | +#!/usr/bin/env python3 |
| 2 | +""" |
| 3 | +Generate 1M peptide datasets using ESM2 for algorithm benchmarking. |
| 4 | +Creates 4 datasets using ESM2 masked language model: 8mer, 9mer, 10mer, 11mer. |
| 5 | +Much faster than ProtGPT-2 and doesn't require fake proteome generation. |
| 6 | +""" |
| 7 | + |
| 8 | +import torch |
| 9 | +import sys |
| 10 | +import time |
| 11 | +import random |
| 12 | +import numpy as np |
| 13 | +from pathlib import Path |
| 14 | +from typing import List, Tuple |
| 15 | +from tqdm import tqdm |
| 16 | +from transformers import EsmTokenizer, EsmForMaskedLM |
| 17 | + |
| 18 | +def setup_esm2_model(model_name: str = "facebook/esm2_t6_8M_UR50D") -> Tuple[EsmTokenizer, EsmForMaskedLM]: |
| 19 | + """ |
| 20 | + Set up ESM2 model and tokenizer. |
| 21 | + Using smallest model (8M) for speed - can upgrade to larger models if needed. |
| 22 | + """ |
| 23 | + print(f"🔄 Loading ESM2 model: {model_name}") |
| 24 | + |
| 25 | + tokenizer = EsmTokenizer.from_pretrained(model_name) |
| 26 | + model = EsmForMaskedLM.from_pretrained(model_name) |
| 27 | + |
| 28 | + # Move to GPU if available |
| 29 | + device = torch.device('cuda' if torch.cuda.is_available() else 'cpu') |
| 30 | + model.to(device) |
| 31 | + model.eval() |
| 32 | + |
| 33 | + print(f"✅ ESM2 model loaded on {device}") |
| 34 | + return tokenizer, model |
| 35 | + |
| 36 | +def generate_peptide_with_esm2(tokenizer, model, length: int, seed_sequence: str = None) -> str: |
| 37 | + """ |
| 38 | + Generate a single peptide of specified length using ESM2. |
| 39 | + Uses iterative masking and prediction. |
| 40 | + """ |
| 41 | + device = next(model.parameters()).device |
| 42 | + |
| 43 | + # Start with a random seed or provided sequence |
| 44 | + if seed_sequence is None: |
| 45 | + # Create random starting sequence |
| 46 | + amino_acids = "ACDEFGHIKLMNPQRSTVWY" |
| 47 | + sequence = ''.join(random.choices(amino_acids, k=length)) |
| 48 | + else: |
| 49 | + sequence = seed_sequence[:length].ljust(length, 'A') # Pad or truncate to desired length |
| 50 | + |
| 51 | + # Iteratively improve the sequence using ESM2 |
| 52 | + num_iterations = 3 # Number of refinement passes |
| 53 | + |
| 54 | + for iteration in range(num_iterations): |
| 55 | + # Randomly mask 1-2 positions |
| 56 | + masked_sequence = list(sequence) |
| 57 | + mask_positions = random.sample(range(length), min(2, length // 3)) |
| 58 | + |
| 59 | + for pos in mask_positions: |
| 60 | + masked_sequence[pos] = tokenizer.mask_token |
| 61 | + |
| 62 | + masked_text = ''.join(masked_sequence) |
| 63 | + |
| 64 | + # Tokenize and predict |
| 65 | + inputs = tokenizer(masked_text, return_tensors="pt").to(device) |
| 66 | + |
| 67 | + with torch.no_grad(): |
| 68 | + outputs = model(**inputs) |
| 69 | + predictions = outputs.logits |
| 70 | + |
| 71 | + # Get predictions for masked positions |
| 72 | + new_sequence = list(sequence) |
| 73 | + for pos in mask_positions: |
| 74 | + token_id = inputs.input_ids[0, pos + 1] # +1 for CLS token |
| 75 | + predicted_token_id = torch.argmax(predictions[0, pos + 1]).item() |
| 76 | + predicted_token = tokenizer.decode([predicted_token_id]) |
| 77 | + |
| 78 | + # Only use valid amino acids |
| 79 | + if predicted_token in "ACDEFGHIKLMNPQRSTVWY": |
| 80 | + new_sequence[pos] = predicted_token |
| 81 | + |
| 82 | + sequence = ''.join(new_sequence) |
| 83 | + |
| 84 | + return sequence |
| 85 | + |
| 86 | +def generate_esm2_peptides(tokenizer, model, length: int, count: int, batch_size: int = 100) -> List[str]: |
| 87 | + """Generate multiple peptides of specified length using ESM2.""" |
| 88 | + peptides = set() # Use set to avoid duplicates |
| 89 | + |
| 90 | + print(f"🔄 Generating {count:,} unique {length}-mer peptides...") |
| 91 | + |
| 92 | + with tqdm(total=count, desc=f"ESM2 {length}mers") as pbar: |
| 93 | + while len(peptides) < count: |
| 94 | + batch_peptides = [] |
| 95 | + |
| 96 | + # Generate batch |
| 97 | + for _ in range(min(batch_size, count - len(peptides))): |
| 98 | + peptide = generate_peptide_with_esm2(tokenizer, model, length) |
| 99 | + batch_peptides.append(peptide) |
| 100 | + |
| 101 | + # Add unique peptides |
| 102 | + initial_size = len(peptides) |
| 103 | + peptides.update(batch_peptides) |
| 104 | + new_peptides = len(peptides) - initial_size |
| 105 | + |
| 106 | + pbar.update(new_peptides) |
| 107 | + |
| 108 | + return list(peptides)[:count] |
| 109 | + |
| 110 | +def write_fasta(sequences: List[str], output_file: Path, prefix: str = "peptide"): |
| 111 | + """Write peptide sequences to FASTA format.""" |
| 112 | + with open(output_file, 'w') as f: |
| 113 | + for i, seq in enumerate(sequences, 1): |
| 114 | + f.write(f">{prefix}_{i:07d}\n{seq}\n") |
| 115 | + |
| 116 | +def generate_esm2_datasets(): |
| 117 | + """Generate all 4 ESM2-based peptide datasets.""" |
| 118 | + |
| 119 | + # Configuration |
| 120 | + base_dir = Path("/Users/chris/Desktop/Griffith Lab/Peptide Sequence Synthesis") |
| 121 | + output_dir = base_dir / "data" / "ESM2_1M_Peptides" |
| 122 | + |
| 123 | + lengths = [8, 9, 10, 11] |
| 124 | + count = 1_000_000 # 1 million peptides per dataset |
| 125 | + |
| 126 | + # Ensure output directory exists |
| 127 | + output_dir.mkdir(parents=True, exist_ok=True) |
| 128 | + |
| 129 | + print(f"🧬 Generating ESM2-based 1M peptide datasets") |
| 130 | + print(f"📁 Output directory: {output_dir}") |
| 131 | + print(f"🔢 Lengths: {lengths}") |
| 132 | + print(f"📈 Count per dataset: {count:,}") |
| 133 | + print("=" * 60) |
| 134 | + |
| 135 | + # Set up ESM2 model |
| 136 | + tokenizer, model = setup_esm2_model() |
| 137 | + |
| 138 | + total_start_time = time.time() |
| 139 | + |
| 140 | + # Generate datasets for each length |
| 141 | + for length in lengths: |
| 142 | + print(f"\n🧪 Generating ESM2 {length}-mer dataset...") |
| 143 | + start_time = time.time() |
| 144 | + |
| 145 | + peptides = generate_esm2_peptides(tokenizer, model, length, count) |
| 146 | + |
| 147 | + # Save to file |
| 148 | + output_file = output_dir / f"esm2_{length}mer_1M.fasta" |
| 149 | + write_fasta(peptides, output_file, prefix=f"esm2_{length}mer") |
| 150 | + |
| 151 | + elapsed = time.time() - start_time |
| 152 | + print(f"✅ Saved {len(peptides):,} unique {length}-mer peptides to {output_file.name}") |
| 153 | + print(f"⏱️ Time: {elapsed:.1f} seconds ({elapsed/60:.1f} minutes)") |
| 154 | + |
| 155 | + # Optional: Save memory by clearing cache |
| 156 | + if torch.cuda.is_available(): |
| 157 | + torch.cuda.empty_cache() |
| 158 | + |
| 159 | + total_elapsed = time.time() - total_start_time |
| 160 | + |
| 161 | + print("\n" + "=" * 60) |
| 162 | + print(f"🎉 ALL ESM2 DATASETS GENERATED SUCCESSFULLY!") |
| 163 | + print(f"⏱️ Total time: {total_elapsed:.1f} seconds ({total_elapsed/60:.1f} minutes)") |
| 164 | + print(f"📁 Output directory: {output_dir}") |
| 165 | + print("\n📋 Generated files:") |
| 166 | + |
| 167 | + # List all generated files |
| 168 | + for length in lengths: |
| 169 | + esm2_file = output_dir / f"esm2_{length}mer_1M.fasta" |
| 170 | + print(f" • {esm2_file.name}") |
| 171 | + |
| 172 | +if __name__ == "__main__": |
| 173 | + # Set random seeds for reproducibility |
| 174 | + random.seed(42) |
| 175 | + np.random.seed(42) |
| 176 | + torch.manual_seed(42) |
| 177 | + |
| 178 | + generate_esm2_datasets() |
0 commit comments