مدیریت توالی‌های DNA و RNA با Biopython.Seq: مفاهیم و کدها

در دوران کنونی که حجم داده‌های بیولوژیکی به طور فزاینده‌ای در حال رشد است، نیاز به ابزارهای محاسباتی قوی برای مدیریت و تحلیل این داده‌ها بیش از پیش احساس می‌شود. در این میان، پایتون به عنوان یک زبان برنامه‌نویسی همه‌کاره و قدرتمند، جایگاه ویژه‌ای در حوزه بیوانفورماتیک پیدا کرده است. Biopython، یک مجموعه کتابخانه‌ای جامع برای پایتون، ابزارهای بی‌نظیری را برای کار با داده‌های بیولوژیکی، از جمله توالی‌های DNA، RNA و پروتئین، فراهم می‌کند. قلب تپنده Biopython برای کار با توالی‌ها، ماژول `Biopython.Seq` است که به کاربران اجازه می‌دهد تا توالی‌ها را به عنوان اشیاء پایتون در نظر بگیرند و عملیات پیچیده‌ای را با سادگی بر روی آن‌ها اعمال کنند.

این مقاله با هدف ارائه یک راهنمای جامع و تخصصی برای مدیریت توالی‌های DNA و RNA با استفاده از `Biopython.Seq` نوشته شده است. ما از مفاهیم پایه‌ای تا تکنیک‌های پیشرفته‌تر را پوشش خواهیم داد، با تاکید بر کاربردها و ارائه مثال‌های کد عملی. مخاطبان این مقاله، متخصصان بیوانفورماتیک، زیست‌شناسان مولکولی، ژنتیک‌دانان و هر محققی هستند که به دنبال راه‌حل‌های برنامه‌نویسی برای تحلیل توالی‌های زیستی می‌باشند. در این مسیر، با اشیاء کلیدی مانند `Seq`، `MutableSeq` و `SeqRecord` آشنا شده و نحوه انجام عملیات حیاتی مانند مکمل‌سازی، ترجمه، ترانسکریپشن و مدیریت فرمت‌های فایلی رایج مانند FASTA و GenBank را خواهیم آموخت. هدف نهایی، توانمندسازی شما برای مدیریت و آنالیز دقیق‌تر و کارآمدتر داده‌های توالی با بهره‌گیری از قدرت بی‌پایان پایتون و Biopython است.

مبانی Biopython.Seq: شیء Seq و انواع توالی

در هسته ماژول `Biopython.Seq`، شیء `Seq` قرار دارد که نمایش‌دهنده یک توالی زیستی (DNA، RNA یا پروتئین) است. این شیء یک کلاس اختصاصی برای نگهداری و دستکاری توالی‌ها فراهم می‌کند که نه تنها رشته نوکلئوتیدی یا آمینواسیدی را ذخیره می‌کند، بلکه توانایی انجام عملیات بیولوژیکی خاص را نیز داراست. درک شیء `Seq` اولین گام برای کار موثر با Biopython در تحلیل توالی است.

ساخت و معرفی شیء Seq

ساده‌ترین راه برای ساخت یک شیء `Seq`، انتقال یک رشته پایتون به سازنده `Seq` است. Biopython به طور خودکار نوع توالی را (DNA، RNA، پروتئین) بر اساس کاراکترهای موجود در رشته حدس نمی‌زند، اما می‌توانیم به صراحت آن را با استفاده از یک الفبای (Alphabet) مشخص کنیم. هرچند، در نسخه‌های جدید Biopython، مفهوم `Alphabet` تا حد زیادی منسوخ شده و `Seq` به طور پیش‌فرض رفتاری عمومی دارد مگر اینکه عملیات خاص بیولوژیکی مانند `translate` یا `complement` فراخوانی شود که در آن زمان ماهیت توالی اهمیت پیدا می‌کند.

یک شیء `Seq` از نوع “immutable” است، به این معنی که پس از ایجاد، نمی‌توان محتوای آن را تغییر داد. این ویژگی ایمنی داده‌ها را تضمین می‌کند و از خطاهای ناخواسته جلوگیری می‌نماید. اگر نیاز به تغییر توالی باشد، باید از شیء `MutableSeq` استفاده کرد که در ادامه به آن خواهیم پرداخت.

from Bio.Seq import Seq

# ایجاد یک توالی DNA
dna_seq = Seq("ATGCGTACGTACGTACGTACGTACGTACG")
print(f"توالی DNA: {dna_seq}")
print(f"نوع شیء: {type(dna_seq)}")

# ایجاد یک توالی RNA (معمولاً با تبدیل DNA به RNA)
rna_seq = Seq("AUGCGUACGACGACGUACGUACGUACG")
print(f"توالی RNA: {rna_seq}")

# ایجاد یک توالی پروتئین
protein_seq = Seq("MCMADRTNRSK")
print(f"توالی پروتئین: {protein_seq}")

تفاوت‌ها و کاربردها: Seq و MutableSeq

همانطور که اشاره شد، شیء `Seq` از نوع “immutable” است. این بدان معناست که اگر بخواهیم یک نوکلئوتید را تغییر دهیم یا بخشی از توالی را حذف کنیم، نمی‌توانیم این کار را مستقیماً روی شیء `Seq` انجام دهیم. برای سناریوهایی که نیاز به تغییر توالی وجود دارد، Biopython شیء `MutableSeq` را ارائه می‌دهد.

`MutableSeq` همانند یک لیست پایتون قابل تغییر است، اما همچنان خصوصیات یک توالی بیولوژیکی را حفظ می‌کند. این ویژگی برای عملیاتی مانند جهش‌زایی (mutagenesis) در silico یا ویرایش توالی‌ها بسیار مفید است.

from Bio.Seq import Seq
from Bio.Seq import MutableSeq

# توالی DNA اولیه (Seq)
dna_seq = Seq("ATGC")
# dna_seq[0] = 'T'  # این خط منجر به خطا می‌شود (TypeError)

# ایجاد یک توالی قابل تغییر (MutableSeq)
mutable_dna_seq = MutableSeq("ATGC")
print(f"توالی قابل تغییر اولیه: {mutable_dna_seq}")

# تغییر یک نوکلئوتید
mutable_dna_seq[0] = 'T'
print(f"توالی پس از تغییر نوکلئوتید اول: {mutable_dna_seq}")

# افزودن نوکلئوتید
mutable_dna_seq.append('A')
print(f"توالی پس از افزودن: {mutable_dna_seq}")

# حذف نوکلئوتید
del mutable_dna_seq[2]
print(f"توالی پس از حذف: {mutable_dna_seq}")

# تبدیل MutableSeq به Seq
final_seq = Seq(mutable_dna_seq)
print(f"توالی نهایی (Seq): {final_seq}")

انتخاب بین `Seq` و `MutableSeq` بستگی به نیاز عملیاتی شما دارد. برای اکثر تحلیل‌های غیرتغییری، `Seq` کافی است و ایمنی بیشتری را فراهم می‌کند. در صورت نیاز به ویرایش، `MutableSeq` ابزار مناسبی است.

عملیات اساسی بر روی توالی‌ها: دستکاری و تحلیل اولیه

پس از درک چگونگی ساخت توالی‌ها، قدم بعدی آشنایی با عملیات اساسی است که می‌توانیم بر روی آن‌ها انجام دهیم. این عملیات شامل دسترسی به عناصر، برش توالی، الحاق، محاسبه طول و جستجوی زیرتوالی‌ها می‌شود که ستون فقرات هر تحلیل بیوانفورماتیکی را تشکیل می‌دهند.

دسترسی به عناصر توالی (Indexing و Slicing)

اشیای `Seq` در Biopython رفتاری مشابه رشته‌های پایتون دارند. این بدان معناست که می‌توانیم با استفاده از اندیس‌گذاری (indexing) به نوکلئوتیدها یا اسیدهای آمینه منفرد دسترسی پیدا کنیم و با استفاده از برش‌دهی (slicing) زیرتوالی‌ها را استخراج کنیم. اندیس‌گذاری در پایتون از صفر شروع می‌شود.

from Bio.Seq import Seq

dna_seq = Seq("AGCTAGCTAGCT")

# دسترسی به نوکلئوتید در یک موقعیت خاص (اندیس 0)
first_base = dna_seq[0]
print(f"نوکلئوتید اول: {first_base}")

# دسترسی به نوکلئوتید در یک موقعیت دیگر (اندیس 5)
sixth_base = dna_seq[5]
print(f"نوکلئوتید ششم: {sixth_base}")

# برش‌دهی برای استخراج زیرتوالی (اندیس 0 تا 4، به اضافه نوکلئوتید در اندیس 4)
sub_seq_1 = dna_seq[0:5]
print(f"زیرتوالی از اندیس 0 تا 4: {sub_seq_1}")

# برش‌دهی از یک نقطه تا انتها
sub_seq_2 = dna_seq[6:]
print(f"زیرتوالی از اندیس 6 تا انتها: {sub_seq_2}")

# برش‌دهی از ابتدا تا یک نقطه
sub_seq_3 = dna_seq[:5]
print(f"زیرتوالی از ابتدا تا اندیس 4: {sub_seq_3}")

# برش‌دهی با گام (هر دو نوکلئوتید)
stepped_seq = dna_seq[::2]
print(f"زیرتوالی با گام 2: {stepped_seq}")

# برش‌دهی معکوس توالی
reversed_seq = dna_seq[::-1]
print(f"توالی معکوس شده: {reversed_seq}")

الحاق توالی‌ها (Concatenation)

الحاق یا چسباندن دو یا چند توالی به یکدیگر، یک عملیات رایج است. اشیای `Seq` را می‌توان با استفاده از عملگر `+` مانند رشته‌های پایتون به هم چسباند. این عملیات یک شیء `Seq` جدید ایجاد می‌کند.

from Bio.Seq import Seq

seq1 = Seq("ATGC")
seq2 = Seq("GTAC")
seq3 = Seq("CCTA")

# الحاق دو توالی
combined_seq = seq1 + seq2
print(f"توالی‌های چسبیده: {combined_seq}")

# الحاق چند توالی
all_combined_seq = seq1 + seq2 + seq3
print(f"چند توالی چسبیده: {all_combined_seq}")

# الحاق با رشته پایتون
new_combined_seq = seq1 + "TTTT"
print(f"الحاق با رشته پایتون: {new_combined_seq}")

طول توالی

برای به دست آوردن طول یک توالی (تعداد نوکلئوتیدها یا اسیدهای آمینه)، از تابع داخلی `len()` پایتون استفاده می‌شود که روی اشیای `Seq` به خوبی کار می‌کند.

from Bio.Seq import Seq

dna_seq = Seq("ATGCGTACGTACGTACGTACGTACGTACG")
rna_seq = Seq("AUGCAUGC")
protein_seq = Seq("MKLGDF")

print(f"طول توالی DNA: {len(dna_seq)}")
print(f"طول توالی RNA: {len(rna_seq)}")
print(f"طول توالی پروتئین: {len(protein_seq)}")

شمارش نوکلئوتیدها/اسیدهای آمینه

شیء `Seq` متد `count()` را برای شمارش تعداد دفعات رخداد یک زیرتوالی یا کاراکتر خاص فراهم می‌کند. این متد بسیار مفید است، به عنوان مثال، برای محاسبه محتوای GC در توالی‌های DNA.

from Bio.Seq import Seq

dna_seq = Seq("ATGCGTACGTACGTACGTACGTACGTACG")

# شمارش تعداد 'G'
g_count = dna_seq.count('G')
print(f"تعداد G در توالی DNA: {g_count}")

# شمارش تعداد 'C'
c_count = dna_seq.count('C')
print(f"تعداد C در توالی DNA: {c_count}")

# محاسبه محتوای GC
gc_content = (g_count + c_count) / len(dna_seq) * 100
print(f"محتوای GC: {gc_content:.2f}%")

# شمارش تعداد یک زیرتوالی
sub_seq_count = dna_seq.count('TAC')
print(f"تعداد زیرتوالی 'TAC': {sub_seq_count}")

جستجوی زیرتوالی‌ها

برای یافتن موقعیت یک زیرتوالی در توالی اصلی، می‌توان از متدهای `find()` و `rfind()` (برای یافتن اولین و آخرین رخداد) یا `index()` و `rindex()` استفاده کرد. تفاوت اصلی بین `find()` و `index()` این است که `find()` در صورت عدم یافتن زیرتوالی، -1 را برمی‌گرداند، در حالی که `index()` یک `ValueError` ایجاد می‌کند.

from Bio.Seq import Seq

dna_seq = Seq("ATGCGTACGTACGTACGTACGTACGTACG")

# یافتن اولین رخداد 'TAC'
first_occurrence = dna_seq.find('TAC')
print(f"اولین رخداد 'TAC' در اندیس: {first_occurrence}")

# یافتن آخرین رخداد 'TAC'
last_occurrence = dna_seq.rfind('TAC')
print(f"آخرین رخداد 'TAC' در اندیس: {last_occurrence}")

# تلاش برای یافتن زیرتوالی ناموجود با find
not_found_find = dna_seq.find('XXX')
print(f"جستجوی 'XXX' با find: {not_found_find}")

# تلاش برای یافتن زیرتوالی ناموجود با index (خطا ایجاد می‌کند)
try:
    not_found_index = dna_seq.index('XXX')
    print(f"جستجوی 'XXX' با index: {not_found_index}")
except ValueError as e:
    print(f"خطا در جستجوی 'XXX' با index: {e}")

با تسلط بر این عملیات اساسی، شما آماده‌اید تا به سراغ توابع بیولوژیکی پیچیده‌تری بروید که Biopython ارائه می‌دهد.

توابع کلیدی برای توالی‌های DNA و RNA: تحولات بیولوژیکی

Biopython فراتر از دستکاری رشته‌ای ساده، عملیات بیولوژیکی اساسی را که برای توالی‌های DNA و RNA حیاتی هستند، مستقیماً بر روی شیء `Seq` فراهم می‌کند. این توابع شامل مکمل‌سازی، معکوس‌مکمل‌سازی، ترانسکریپشن، معکوس ترانسکریپشن و ترجمه می‌شوند که هر یک کاربردهای گسترده‌ای در زیست‌شناسی مولکولی و ژنتیک دارند.

مکمل‌سازی (Complement)

مکمل یک توالی DNA به معنای جایگزینی هر باز با باز جفت‌شونده خود است (A با T، T با A، G با C، C با G). این عملیات در بیولوژی برای درک ساختار دو رشته‌ای DNA و در تکنیک‌هایی مانند PCR و هیبریداسیون اهمیت دارد.

Biopython متد `complement()` را برای این منظور ارائه می‌دهد. این متد یک شیء `Seq` جدید را با توالی مکمل برمی‌گرداند.

from Bio.Seq import Seq

dna_seq = Seq("ATGCGTACGT")
complement_seq = dna_seq.complement()
print(f"توالی DNA اصلی: {dna_seq}")
print(f"توالی مکمل: {complement_seq}")

مکمل معکوس (Reverse Complement)

مکمل معکوس یک توالی، توالی مکملی است که سپس به صورت معکوس خوانده می‌شود (از انتها به ابتدا). این عملیات برای کارهای مرتبط با رشته مکمل DNA در جهت ‘3 به ‘5 و به خصوص در طراحی پرایمرها (primers) برای PCR، نقشه‌برداری ژنوم و سنتز الیگونوکلئوتیدها بسیار حیاتی است.

متد `reverse_complement()` این عمل را انجام می‌دهد:

from Bio.Seq import Seq

dna_seq = Seq("ATGCGTACGT")
reverse_complement_seq = dna_seq.reverse_complement()
print(f"توالی DNA اصلی: {dna_seq}")
print(f"توالی مکمل معکوس: {reverse_complement_seq}")

ترانسکریپشن (Transcription): DNA به RNA

ترانسکریپشن فرآیند رونویسی اطلاعات ژنتیکی از یک رشته DNA به یک رشته mRNA (مولکول RNA پیام‌رسان) است. در این فرآیند، تیمین (T) در DNA با یوراسیل (U) در RNA جایگزین می‌شود و بقیه بازها (A, G, C) ثابت می‌مانند. Biopython با متد `transcribe()` این فرآیند را مدل‌سازی می‌کند.

from Bio.Seq import Seq

dna_template = Seq("ATGCGTACGTACGTACGTACGTACGTACG")
mrna_seq = dna_template.transcribe()
print(f"رشته DNA الگو: {dna_template}")
print(f"رشته mRNA رونویسی شده: {mrna_seq}")

معکوس ترانسکریپشن (Back-transcription): RNA به DNA

معکوس ترانسکریپشن، فرآیند برعکس ترانسکریپشن است که در آن یک رشته DNA از یک الگوی RNA سنتز می‌شود. این فرآیند در ویروس‌های دارای ژنوم RNA (مانند HIV) و در تکنیک‌های آزمایشگاهی مانند RT-PCR (Reverse Transcription PCR) مهم است. Biopython متد `back_transcribe()` را برای این منظور فراهم می‌کند.

from Bio.Seq import Seq

rna_seq = Seq("AUGCAUGCGUACGUACGUACGUACGUACG")
cdna_seq = rna_seq.back_transcribe()
print(f"رشته RNA: {rna_seq}")
print(f"رشته cDNA (معکوس ترانسکریپشن شده): {cdna_seq}")

ترجمه (Translation): RNA به پروتئین

ترجمه فرآیندی است که در آن اطلاعات ژنتیکی موجود در mRNA به یک زنجیره پلی‌پپتیدی (پروتئین) تبدیل می‌شود. این فرآیند از طریق کدون‌ها (سه‌تایی‌های نوکلئوتیدی) و با استفاده از جدول کد ژنتیکی انجام می‌شود. Biopython متد `translate()` را برای انجام این عملیات ارائه می‌دهد.

متد `translate()` به صورت پیش‌فرض از جدول کد ژنتیکی استاندارد استفاده می‌کند، اما می‌توان جدول‌های کد ژنتیکی جایگزین را نیز مشخص کرد (مثلاً برای میتوکندری یا باکتری‌ها) با استفاده از پارامتر `table`.

from Bio.Seq import Seq
from Bio.Data import CodonTable

# توالی mRNA
mrna_seq = Seq("AUGCAUGCGUACGUACGUACGUACGUACGUGAA")
protein_seq = mrna_seq.translate()
print(f"رشته mRNA: {mrna_seq}")
print(f"رشته پروتئین ترجمه شده (استاندارد): {protein_seq}")

# ترجمه با در نظر گرفتن کدون‌های پایان (stop codons) به عنوان X
protein_with_stops = mrna_seq.translate(to_stop=False)
print(f"رشته پروتئین ترجمه شده (با stop به عنوان X): {protein_with_stops}")

# ترجمه با جدول کد ژنتیکی جایگزین (مثلاً میتوکندری مهره‌داران)
# برای مشاهده لیست جدول‌ها: dir(CodonTable.unambiguous_dna_by_name)
# یا CodonTable.ambiguous_dna_by_name
# جدول 2 (Vertebrate Mitochondrial)
mitochondrial_protein_seq = mrna_seq.translate(table=2)
print(f"رشته پروتئین ترجمه شده (میتوکندری مهره‌داران): {mitochondrial_protein_seq}")

# مثال دیگر برای translation
dna_seq_for_translation = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG")
# ابتدا DNA را به RNA تبدیل می‌کنیم
mrna_from_dna = dna_seq_for_translation.transcribe()
translated_protein = mrna_from_dna.translate()
print(f"\nرشته DNA برای ترجمه: {dna_seq_for_translation}")
print(f"رشته mRNA حاصل: {mrna_from_dna}")
print(f"رشته پروتئین ترجمه شده: {translated_protein}")

# بررسی ترجمه فریم‌های مختلف
# فریم 0
print(f"فریم 0: {dna_seq_for_translation[0:].translate(to_stop=True)}")
# فریم 1
print(f"فریم 1: {dna_seq_for_translation[1:].translate(to_stop=True)}")
# فریم 2
print(f"فریم 2: {dna_seq_for_translation[2:].translate(to_stop=True)}")

پیدا کردن ORF (Open Reading Frame)

Open Reading Frame (ORF) یا قاب خوانش باز، یک توالی پیوسته از کدون‌ها در یک مولکول DNA یا RNA است که با یک کدون آغازگر (معمولاً AUG) شروع شده و با یک کدون پایان (UAA, UAG, UGA) خاتمه می‌یابد و پتانسیل کد کردن یک پروتئین را دارد. پیدا کردن ORFها یک گام حیاتی در پیش‌بینی ژن‌ها و پروتئین‌ها است.

Biopython مستقیماً یک متد `find_orfs()` بر روی شیء `Seq` ندارد، اما با ترکیب `translate()` و برش‌دهی (slicing) می‌توانیم ORFها را جستجو کنیم. معمولاً برای این کار باید هر سه قاب خوانش (reading frame) در هر دو رشته (حس و پادحس) را بررسی کنیم. این یک مثال ساده است و ابزارهای پیشرفته‌تر (مانند EMBOSS getorf یا نرم‌افزارهای دسکتاپ) برای شناسایی دقیق ORFها در ژنوم‌های بزرگ وجود دارند.

from Bio.Seq import Seq

# یک توالی DNA فرضی
dna_seq_with_orf = Seq("ATGCGGATGCAGCATGCGTGAAGCAATAGCCGATGCAATAATAA")

# ترجمه در سه فریم خوانش (reading frames)
print(f"توالی DNA: {dna_seq_with_orf}")

# فریم 0 (از اندیس 0 شروع می‌شود)
protein_frame0 = dna_seq_with_orf.translate(table=1, to_stop=True)
print(f"فریم 0 پروتئین: {protein_frame0}")

# فریم 1 (از اندیس 1 شروع می‌شود)
protein_frame1 = dna_seq_with_orf[1:].translate(table=1, to_stop=True)
print(f"فریم 1 پروتئین: {protein_frame1}")

# فریم 2 (از اندیس 2 شروع می‌شود)
protein_frame2 = dna_seq_with_orf[2:].translate(table=1, to_stop=True)
print(f"فریم 2 پروتئین: {protein_frame2}")

# بررسی رشته مکمل معکوس نیز برای ORF ها
reverse_complement_dna = dna_seq_with_orf.reverse_complement()
print(f"\nتوالی مکمل معکوس: {reverse_complement_dna}")

protein_frame0_rc = reverse_complement_dna.translate(table=1, to_stop=True)
print(f"فریم 0 پروتئین (مکمل معکوس): {protein_frame0_rc}")

protein_frame1_rc = reverse_complement_dna[1:].translate(table=1, to_stop=True)
print(f"فریم 1 پروتئین (مکمل معکوس): {protein_frame1_rc}")

protein_frame2_rc = reverse_complement_dna[2:].translate(table=1, to_stop=True)
print(f"فریم 2 پروتئین (مکمل معکوس): {protein_frame2_rc}")

# تابع ساده برای یافتن ORFهای بالای طول مشخص
def find_all_orfs(sequence, min_protein_length=10):
    orfs = []
    # برای هر دو رشته (حس و پادحس)
    for strand, dna in [(+1, sequence), (-1, sequence.reverse_complement())]:
        for frame in range(3):
            protein = dna[frame:].translate(table=1, to_stop=False)
            current_orf = ""
            start_pos = frame
            for i, aa in enumerate(protein):
                if aa == 'M' and not current_orf: # شروع ORF
                    current_orf = str(protein[i:])
                    start_pos_in_dna = frame + i * 3
                if aa == '*' and current_orf: # پایان ORF
                    end_pos_in_dna = start_pos_in_dna + len(current_orf[:-1]) * 3 + 3
                    if len(current_orf[:-1]) >= min_protein_length:
                        orfs.append({
                            'sequence': current_orf[:-1],
                            'strand': strand,
                            'frame': frame,
                            'start_codon_pos': start_pos_in_dna,
                            'end_codon_pos': end_pos_in_dna
                        })
                    current_orf = ""
            # اگر ORF به پایان نرسیده و طول مناسب دارد
            if current_orf and len(current_orf[:-1]) >= min_protein_length:
                orfs.append({
                    'sequence': current_orf[:-1],
                    'strand': strand,
                    'frame': frame,
                    'start_codon_pos': start_pos_in_dna,
                    'end_codon_pos': start_pos_in_dna + len(current_orf) * 3
                })
    return orfs

orfs_found = find_all_orfs(dna_seq_with_orf, min_protein_length=4)
print("\nORF های یافت شده (طول پروتئین حداقل 4):")
for orf in orfs_found:
    print(f"  رشته: {orf['sequence']}, جهت: {orf['strand']}, فریم: {orf['frame']}, "
          f"شروع: {orf['start_codon_pos']}, پایان: {orf['end_codon_pos']}")

توابع بالا، سنگ بنای بسیاری از تحلیل‌های بیولوژیکی هستند و Biopython آن‌ها را به سادگی در دسترس قرار می‌دهد.

کار با توالی‌ها و Annotation ها: معرفی SeqRecord و SeqFeature

توالی‌های DNA، RNA یا پروتئین تنها رشته‌ای از کاراکترها نیستند؛ آن‌ها اغلب با اطلاعات اضافی ارزشمندی همراه هستند. این اطلاعات می‌توانند شامل شناسه توالی، نام، توضیحات، مبدا، و همچنین ویژگی‌های ژنتیکی مانند جایگاه ژن‌ها، مناطق کدشونده (CDS)، نواحی تنظیمی و سایر اطلاعات ساختاری باشند. شیء `Seq` در Biopython به تنهایی برای نگهداری این ابرداده (metadata) کافی نیست.

برای پاسخگویی به این نیاز، Biopython اشیاء `SeqRecord` و `SeqFeature` را معرفی می‌کند. `SeqRecord` به عنوان یک “رکورد توالی” عمل می‌کند که یک شیء `Seq` را به همراه تمام ابرداده‌های مرتبط خود دربرمی‌گیرد. `SeqFeature` نیز برای حاشیه‌نویسی (annotation) مناطق خاصی از توالی استفاده می‌شود.

چرا فقط Seq کافی نیست؟ نیاز به Metadata

فرض کنید فایلی با فرمت FASTA را می‌خوانید. هر رکورد در یک فایل FASTA علاوه بر توالی، یک خط سربرگ (header) نیز دارد که شامل شناسه، نام و توضیحاتی در مورد توالی است. شیء `Seq` تنها خود توالی را نگه می‌دارد و راهی برای ذخیره این اطلاعات اضافی ندارد. `SeqRecord` دقیقا برای همین منظور طراحی شده است.

from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord

# ایجاد یک توالی ساده
my_seq = Seq("ATGCGTACGTACGTACGTACG")

# ایجاد یک SeqRecord از توالی
# شامل: id, name, description
my_record = SeqRecord(
    my_seq,
    id="NC_000913.3",
    name="E.coli K-12 MG1655",
    description="Complete genome"
)

print(f"شناسه (ID): {my_record.id}")
print(f"نام (Name): {my_record.name}")
print(f"توضیحات (Description): {my_record.description}")
print(f"توالی (Seq): {my_record.seq}")
print(f"طول توالی: {len(my_record.seq)}")

افزودن Annotation ها و Features به SeqRecord

`SeqRecord` دارای دو دیکشنری مهم برای ذخیره اطلاعات بیشتر است: `annotations` و `features`.
`annotations` برای اطلاعات کلی در مورد کل توالی مانند منبع (source)، ارگانیسم، مرجع‌های PubMed و غیره استفاده می‌شود.
`features` فهرستی از اشیاء `SeqFeature` است که نشان‌دهنده نواحی مشخصی در توالی مانند ژن‌ها، CDSها، tRNAها و غیره هستند.

from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
from Bio.SeqFeature import SeqFeature, FeatureLocation

my_seq = Seq("ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG")

my_record = SeqRecord(
    my_seq,
    id="gene_example",
    name="Gene_Alpha",
    description="An example gene sequence"
)

# افزودن Annotation ها
my_record.annotations["organism"] = "Homo sapiens"
my_record.annotations["molecule_type"] = "DNA"
my_record.annotations["date"] = "2023-10-27"
my_record.annotations["source"] = "Synthetically derived"

print(f"\nAnnotation ها: {my_record.annotations}")

# افزودن Feature ها
# یک FeatureLocation مشخص کننده محدوده یک Feature است.
# اینجا یک منطقه کدشونده (CDS) از اندیس 0 تا 39 (با رعایت zero-based indexing)
gene_feature = SeqFeature(FeatureLocation(0, 39), type="CDS", id="CDS_1", qualifiers={"product": ["Hypothetical Protein"], "codon_start": [1]})
my_record.features.append(gene_feature)

# افزودن یک Feature دیگر برای یک ناحیه تنظیمی
promoter_feature = SeqFeature(FeatureLocation(5, 15), type="promoter", qualifiers={"note": ["TATA box region"]})
my_record.features.append(promoter_feature)

print(f"\nتعداد Feature ها: {len(my_record.features)}")
for feature in my_record.features:
    print(f"  نوع Feature: {feature.type}")
    print(f"  موقعیت: {feature.location}")
    print(f"  توالی Feature: {feature.extract(my_seq)}")
    print(f"  Qualifier ها: {feature.qualifiers}")

`FeatureLocation` علاوه بر شروع و پایان، می‌تواند اطلاعات مربوط به جهت‌گیری (strand) و موقعیت‌های غیرپیوسته (مثلاً برای اینترون‌ها و اگزون‌ها) را نیز ذخیره کند. این انعطاف‌پذیری `SeqRecord` را به یک ابزار قدرتمند برای کار با داده‌های حاشیه‌نویسی شده مانند فایل‌های GenBank تبدیل می‌کند.

با استفاده از `SeqRecord` و `SeqFeature`، می‌توانید نه تنها توالی‌ها را مدیریت کنید، بلکه تمام اطلاعات بیولوژیکی و متا دیتای مرتبط با آن‌ها را نیز به شکلی ساختاریافته حفظ و پردازش کنید.

ورودی و خروجی توالی‌ها: کار با فایل‌های FASTA و GenBank

در بیوانفورماتیک، توالی‌ها اغلب در قالب فایل‌هایی با فرمت‌های استاندارد ذخیره می‌شوند. دو تا از رایج‌ترین فرمت‌ها، FASTA و GenBank هستند. Biopython ماژول `Bio.SeqIO` را برای خواندن و نوشتن توالی‌ها به این فرمت‌ها و بسیاری فرمت‌های دیگر ارائه می‌دهد. `Bio.SeqIO` ابزاری بسیار قدرتمند و انعطاف‌پذیر برای مدیریت داده‌های توالی از/به فایل‌ها است.

خواندن فایل‌های FASTA با SeqIO.parse

فرمت FASTA یکی از ساده‌ترین و رایج‌ترین فرمت‌ها برای ذخیره توالی‌های بیولوژیکی است. هر رکورد FASTA با یک خط سربرگ شروع می‌شود که با `>` مشخص شده و به دنبال آن توالی در یک یا چند خط قرار می‌گیرد. `Bio.SeqIO.parse()` تابعی است که یک ژنراتور (generator) از اشیاء `SeqRecord` را از یک فایل FASTA بازگردانده و به شما امکان می‌دهد تا توالی‌های موجود در فایل را به صورت کارآمد، حتی در فایل‌های بسیار بزرگ، پردازش کنید.

from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord

# ایجاد یک فایل FASTA نمونه
fasta_content = """
>seq1 | Description for sequence 1
ATGCGTACGTACGTACGTACGTACGTACG
>seq2 | Another sequence example
GATTACA
>seq3 | Yet another sequence
AAAGGGTTTCCC
"""

with open("example.fasta", "w") as f:
    f.write(fasta_content)

# خواندن فایل FASTA
print("خواندن فایل FASTA:")
for record in SeqIO.parse("example.fasta", "fasta"):
    print(f"  ID: {record.id}")
    print(f"  Description: {record.description}")
    print(f"  Sequence: {record.seq}")
    print(f"  طول: {len(record.seq)}")
    print("-" * 20)

# مثال: استخراج توالی‌های با طول خاص
print("\nتوالی‌های با طول بیشتر از 10:")
long_sequences = []
for record in SeqIO.parse("example.fasta", "fasta"):
    if len(record.seq) > 10:
        long_sequences.append(record)

for record in long_sequences:
    print(f"  ID: {record.id}, طول: {len(record.seq)}")

نوشتن توالی‌ها به فرمت FASTA

`SeqIO.write()` به شما امکان می‌دهد تا لیستی از اشیاء `SeqRecord` را به یک فایل با فرمت مشخص (مانند FASTA) بنویسید.

from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord

# ایجاد لیستی از SeqRecord ها
record1 = SeqRecord(
    Seq("ATGCGTACGTACG"),
    id="my_seq_1",
    name="Sequence_A",
    description="A test sequence for writing"
)
record2 = SeqRecord(
    Seq("GATTACA"),
    id="my_seq_2",
    name="Sequence_B",
    description="Another test sequence"
)
records_to_write = [record1, record2]

# نوشتن به فایل FASTA
output_filename = "output.fasta"
count = SeqIO.write(records_to_write, output_filename, "fasta")
print(f"\n{count} رکورد به فایل '{output_filename}' نوشته شد.")

# محتوای فایل output.fasta را چاپ می‌کنیم تا مطمئن شویم
with open(output_filename, "r") as f:
    print("\nمحتوای فایل output.fasta:")
    print(f.read())

خواندن فایل‌های GenBank

فرمت GenBank بسیار غنی‌تر از FASTA است و علاوه بر توالی، حاوی مقدار زیادی ابرداده و حاشیه‌نویسی (annotations و features) است. این شامل اطلاعات مربوط به مبدا ارگانیسم، مقالات مرجع، ژن‌ها، مناطق کدشونده، پروتئین‌های کدگذاری شده و غیره می‌شود. `SeqIO.parse()` می‌تواند فایل‌های GenBank را نیز بخواند و هر رکورد را به یک شیء `SeqRecord` با تمام ابرداده‌های مرتبط تبدیل کند.

from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
from Bio.SeqFeature import SeqFeature, FeatureLocation

# ایجاد یک فایل GenBank نمونه (به دلیل پیچیدگی ساخت دستی، از یک رشته ثابت استفاده می‌کنیم)
genbank_content = """
LOCUS       test_locus              39 bp    DNA     linear   INV 27-OCT-2023
DEFINITION  An example gene sequence for Biopython testing.
ACCESSION   test_id
VERSION     test_id.1
KEYWORDS    Biopython; example; DNA.
SOURCE      Homo sapiens
  ORGANISM  Homo sapiens
            Eukaryota; Metazoa; Chordata; Craniata; Vertebrata; Euteleostomi;
            Mammalia; Eutheria; Euarchontoglires; Primates; Hominidae; Homo.
FEATURES             Location/Qualifiers
     source          1..39
                     /organism="Homo sapiens"
                     /mol_type="genomic DNA"
                     /db_xref="taxon:9606"
     CDS             1..39
                     /codon_start=1
                     /product="Hypothetical Protein"
                     /protein_id="XP_123456.1"
                     /translation="MGDVMAALKRVPDS"
ORIGIN
        1 atggccattg taatgggccg ctgaaaggat gcccgatag
//
"""
with open("example.gb", "w") as f:
    f.write(genbank_content)

# خواندن فایل GenBank
print("\nخواندن فایل GenBank:")
for record in SeqIO.parse("example.gb", "genbank"):
    print(f"  ID: {record.id}")
    print(f"  Name: {record.name}")
    print(f"  Description: {record.description}")
    print(f"  Sequence: {record.seq}")
    print(f"  Annotations: {record.annotations}")
    print(f"  تعداد Feature ها: {len(record.features)}")
    for feature in record.features:
        print(f"    نوع Feature: {feature.type}")
        print(f"    موقعیت Feature: {feature.location}")
        if feature.type == "CDS":
            print(f"    محصول پروتئین: {feature.qualifiers.get('product', ['N/A'])[0]}")
            print(f"    ترجمه پروتئین: {feature.qualifiers.get('translation', ['N/A'])[0]}")
    print("-" * 20)

تفاوت‌های FASTA و GenBank

• **FASTA:** فرمتی ساده و فشرده برای ذخیره توالی‌های خام (DNA، RNA، پروتئین). هر رکورد شامل یک خط سربرگ (شناسه و توضیحات) و سپس خود توالی است. مناسب برای ذخیره توالی‌های پرشمار یا زمانی که ابرداده‌های تفصیلی نیاز نیستند.
• **GenBank:** فرمتی بسیار جامع که علاوه بر توالی، شامل مقادیر زیادی ابرداده ساختاریافته و حاشیه‌نویسی (features و annotations) است. این ابرداده شامل اطلاعات بیولوژیکی، انتشارات مرتبط، ویژگی‌های ژنتیکی (ژن‌ها، مناطق کدشونده، پروموترها و غیره) و جزئیات مربوط به ارگانیسم است. مناسب برای ذخیره توالی‌های دارای حاشیه‌نویسی غنی، مانند ژنوم‌های کامل یا مناطق ژنومی خاص.

`Bio.SeqIO` از فرمت‌های بسیار دیگری نیز پشتیبانی می‌کند (مانند EMBL، SwissProt، PDB و غیره). با استفاده از این ماژول، می‌توانید به راحتی داده‌های توالی را بین فرمت‌های مختلف تبدیل کرده و برای تحلیل‌های بعدی آماده کنید.

سناریوهای پیشرفته و کاربردهای عملی: فراتر از مبانی

پس از آشنایی با مبانی `Biopython.Seq` و `Bio.SeqIO`، اکنون زمان آن است که به برخی از سناریوهای پیشرفته‌تر و کاربردهای عملی بپردازیم که قدرت این کتابخانه را در تحلیل‌های بیوانفورماتیکی واقعی نشان می‌دهند. این موارد شامل تحلیل تفاوت‌های توالی، ادغام با ابزارهای دیگر و استفاده در پایپ‌لاین‌های پیچیده‌تر می‌شود.

تحلیل تفاوت‌های توالی (پیدا کردن SNP ها و جهش‌ها)

یکی از کاربردهای کلیدی در ژنتیک و بیولوژی تکاملی، شناسایی و تحلیل تفاوت‌ها بین توالی‌ها است. این تفاوت‌ها می‌توانند شامل جهش‌های نقطه‌ای (Single Nucleotide Polymorphisms – SNPs)، حذف‌ها (deletions) یا درج‌ها (insertions) باشند. Biopython ابزارهایی برای تراز توالی‌ها (Multiple Sequence Alignment – MSA) ارائه می‌دهد که با ترکیب آن‌ها با قابلیت‌های `Seq` می‌توان به این تحلیل‌ها پرداخت.

from Bio.Seq import Seq
from Bio.Align import PairwiseAligner

# دو توالی DNA برای مقایسه
seq_ref = Seq("ATGCGTACGTACGTACGTACGTACGTACG")
seq_mut = Seq("ATGCTTAGGTACGTACGTTCGTACGTACG")

print(f"توالی مرجع: {seq_ref}")
print(f"توالی جهش‌یافته: {seq_mut}")

# مقایسه نوکلئوتید به نوکلئوتید برای یافتن SNP ها (فرض تراز بودن توالی‌ها)
print("\nتفاوت‌های نوکلئوتید به نوکلئوتید (SNP ها):")
for i in range(min(len(seq_ref), len(seq_mut))):
    if seq_ref[i] != seq_mut[i]:
        print(f"  تفاوت در موقعیت {i}: مرجع {seq_ref[i]} -> جهش‌یافته {seq_mut[i]}")

# استفاده از PairwiseAligner برای یافتن تفاوت‌ها (برای توالی‌های ناتراز)
# این یک تراز محلی یا سراسری را انجام می‌دهد و می‌تواند ایندل‌ها را نیز شناسایی کند.
aligner = PairwiseAligner()
aligner.mode = 'global' # یا 'local'
aligner.match_score = 1
aligner.mismatch_score = -1
aligner.open_gap_score = -2
aligner.extend_gap_score = -0.5

alignments = aligner.align(seq_ref, seq_mut)

# نمایش بهترین تراز
print("\nبهترین تراز سراسری:")
for alignment in sorted(alignments, key=lambda x: x.score, reverse=True)[:1]:
    print(alignment)

# استخراج تفاوت‌ها از تراز (مثالی ساده، در عمل نیاز به پردازش دقیق خروجی aligner دارد)
aligned_seq_ref, aligned_seq_mut = str(alignment[0]), str(alignment[1])
print("\nتفاوت‌ها شامل ایندل‌ها از تراز:")
current_pos_ref = 0
current_pos_mut = 0
for i in range(len(aligned_seq_ref)):
    if aligned_seq_ref[i] != aligned_seq_mut[i]:
        if aligned_seq_ref[i] == '-': # Insertion in seq_mut relative to seq_ref
            print(f"  Insertion در توالی جهش‌یافته در موقعیت {current_pos_mut}: {aligned_seq_mut[i]} (بعد از {current_pos_ref} در مرجع)")
        elif aligned_seq_mut[i] == '-': # Deletion in seq_mut relative to seq_ref
            print(f"  Deletion در توالی جهش‌یافته در موقعیت {current_pos_ref}: {aligned_seq_ref[i]} (از توالی مرجع حذف شده)")
        else: # SNP
            print(f"  SNP در موقعیت مرجع {current_pos_ref} (موقعیت جهش‌یافته {current_pos_mut}): {aligned_seq_ref[i]} -> {aligned_seq_mut[i]}")
    if aligned_seq_ref[i] != '-':
        current_pos_ref += 1
    if aligned_seq_mut[i] != '-':
        current_pos_mut += 1

استفاده در پایپ‌لاین‌های بیوانفورماتیکی

Biopython و به خصوص ماژول `Seq`، اغلب به عنوان بلوک‌های سازنده در پایپ‌لاین‌های بیوانفورماتیکی بزرگتر عمل می‌کنند. مثلاً:

• **طراحی پرایمر:** با استفاده از `reverse_complement()` و جستجوی زیرتوالی، می‌توان پرایمرهای PCR را طراحی و بررسی کرد.
• **آماده‌سازی داده برای BLAST:** توالی‌ها را از فایل‌های مختلف خوانده، فیلتر کرده (مثلاً بر اساس طول یا محتوای GC) و به فرمت FASTA برای ارسال به ابزارهای BLAST آماده کرد.
• **تحلیل‌های فیلوژنتیکی:** توالی‌های پروتئینی را از DNA ترجمه کرده، آن‌ها را تراز کرده و درخت‌های فیلوژنتیکی ساخت.
• **پردازش داده‌های توالی نسل جدید (NGS):** هرچند Biopython مستقیماً برای فایل‌های حجیم NGS (مانند BAM یا VCF) بهینه نیست، اما می‌تواند برای پردازش توالی‌های استخراج شده یا ایجاد توالی‌های مرجع استفاده شود.

from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
import random

# سناریو: فیلتر کردن توالی‌ها و آماده‌سازی برای BLAST
# فرض کنید یک فایل ورودی FASTA داریم
input_fasta_content = """
>geneA | Hypothetical protein A
ATGGGGCTAAACTAGACTGAGATTAGCTAGCTAGTAGCTAG
>geneB | Short sequence, not relevant
ATGC
>geneC | Another important gene
ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG
>geneD | Contains an unknown region
ATGCGTACGTNNNNNNNNNCGTACGTACG
"""
with open("input.fasta", "w") as f:
    f.write(input_fasta_content)

# توالی‌های فیلتر شده بر اساس طول و عدم وجود N ها
filtered_records = []
min_length = 20
max_unknown_bases = 0.1 # حداکثر 10 درصد N در توالی

for record in SeqIO.parse("input.fasta", "fasta"):
    num_n = record.seq.count('N')
    # بررسی طول و تعداد N
    if len(record.seq) >= min_length and (num_n / len(record.seq)) <= max_unknown_bases:
        filtered_records.append(record)
        
output_blast_fasta = "filtered_for_blast.fasta"
count_filtered = SeqIO.write(filtered_records, output_blast_fasta, "fasta")
print(f"\n{count_filtered} رکورد فیلتر شده (طول >= {min_length}, N < {max_unknown_bases*100}%) به '{output_blast_fasta}' نوشته شد.")

with open(output_blast_fasta, "r") as f:
    print("\nمحتوای فایل فیلتر شده برای BLAST:")
    print(f.read())

ادغام با ابزارهای دیگر (مثلاً BLAST)

Biopython می‌تواند به عنوان یک رابط برای ابزارهای خارجی مانند BLAST (Basic Local Alignment Search Tool) عمل کند. این امکان به شما می‌دهد که اسکریپت‌های پایتون خود را با ابزارهای خط فرمان ادغام کنید، بدون نیاز به ترک محیط پایتون.

from Bio.Blast import NCBIWWW
from Bio.Blast import NCBIXML
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord

# این مثال نیاز به اتصال به اینترنت دارد و ممکن است زمان‌بر باشد.
# همچنین نیازمند توالی‌های واقعی برای جستجوی معنی‌دار است.

# یک توالی پروتئینی فرضی
# protein_query = Seq("MSKASAGKKKKAA")
# record = SeqRecord(protein_query, id="query_prot", description="Test protein query")

# print(f"\nجستجوی BLASTp برای توالی: {record.id}")

# # انجام جستجوی BLASTp از طریق NCBIWWW (نیاز به اینترنت)
# # result_handle = NCBIWWW.qblast("blastp", "nr", record.format("fasta"))

# # برای جلوگیری از اجرای زنده در تست، از یک فایل نتیجه BLASTP نمونه استفاده می‌کنیم
# # فرض کنید result_handle محتوای یک XML BLAST است
# # در حالت واقعی، result_handle از NCBIWWW.qblast برمی‌گردد

# # نمونه خروجی XML BLAST (فرض شده)
# blast_xml_sample = """
# 
# 
# 
#   blastp
#   BLASTP 2.13.0+
#   Stephen F. Altschul, Thomas L. Madden, Alejandro A. Schäffer, Jinghui Zhang, Zheng Zhang, Webb Miller, and David J. Lipman (1997), "Gapped BLAST and PSI-BLAST: a new generation of protein database search programs", Nucleic Acids Res. 25:3389-3402.
#   nr
#   query_prot
#   Test protein query
#   13
#   
#     
#       BLOSUM62
#       10
#       11
#       1
#       F
#     
#   
# 
# 
#   1
#   query_prot
#   Test protein query
#   13
# 
# 
#   1
#   gb|AAD15697.1|
#   hypothetical protein [Escherichia coli K-12]
#   AAD15697
#   100
#   
#     
#       1
#       26.5682
#       57
#       0.106869
#       1
#       13
#       1
#       13
#       0
#       0
#       7
#       8
#       0
#       13
#       MSKASAGKKKKAA
#       MSKASA-KKRKAA
#       MSKASA KK+KAA
#     
#   
# 
# 
# 
# 
# 
# """
#
# # برای اجرای واقعی: result_handle = NCBIWWW.qblast("blastp", "nr", record.format("fasta"))
# # و سپس: blast_records = NCBIXML.parse(result_handle)
# # اینجا از StringIO برای شبیه‌سازی فایل استفاده می‌کنیم
# from io import StringIO
# blast_records = NCBIXML.parse(StringIO(blast_xml_sample))
#
# for blast_record in blast_records:
#     for alignment in blast_record.alignments:
#         print(f"  ردیف همترازی: {alignment.title}")
#         for hsp in alignment.hsps:
#             print(f"    E-value: {hsp.expect}")
#             print(f"    query: {hsp.query}")
#             print(f"    match: {hsp.match}")
#             print(f"    hit:   {hsp.sbjct}")
#             print("-" * 10)
#
# print("\n(نکته: مثال BLAST بالا برای جلوگیری از بارگذاری شبکه در زمان اجرا، از یک خروجی XML نمونه استفاده کرده است. "
#       "برای اجرای واقعی، خطوط مربوط به NCBIWWW.qblast را از حالت کامنت خارج کنید.)")

این قابلیت‌ها `Biopython.Seq` را به یک ابزار ضروری در جعبه ابزار هر متخصص بیوانفورماتیک تبدیل می‌کند که به دنبال اتوماسیون و کارایی در تحلیل داده‌های توالی است.

مثال‌های عملی برای تحلیل ژنوم

در نهایت، اجازه دهید یک سناریوی کمی جامع‌تر برای تحلیل ژنوم ارائه دهیم که نشان می‌دهد چگونه می‌توان از ابزارهای Biopython در کنار یکدیگر استفاده کرد.

from Bio import SeqIO
from Bio.Seq import Seq
from Bio.SeqRecord import SeqRecord
from Bio.Data import CodonTable

# سناریو: تحلیل یک قطعه ژنومی فرضی
# فرض کنید یک توالی DNA از یک ژنوم داریم و می‌خواهیم:
# 1. محتوای GC آن را محاسبه کنیم.
# 2. فریم‌های خوانش آن را ترجمه کنیم.
# 3. یک ORF خاص را شناسایی کنیم.

genome_segment_dna = Seq("ATGGCTAGCTAGCATGACTAGCTAGATGACTAGCTAGCATGCATGCATGCATGCGCATGCATGCATGCATGCATGCATGCATGCATGCATGTAGCTAGCTAGCTAGCTAGCTAGGGTCTAGCTAGCTAGCTAGCTAGCTGTAGCC")
genome_segment_record = SeqRecord(
    genome_segment_dna,
    id="chr1_segment_1",
    name="Segment 1",
    description="Hypothetical segment from chromosome 1"
)

print(f"قطعه ژنومی: {genome_segment_record.seq}")
print(f"طول قطعه: {len(genome_segment_record.seq)}")

# 1. محاسبه محتوای GC
gc_count = genome_segment_record.seq.count('G') + genome_segment_record.seq.count('C')
total_bases = len(genome_segment_record.seq)
gc_content_percentage = (gc_count / total_bases) * 100
print(f"محتوای GC: {gc_content_percentage:.2f}%")

# 2. ترجمه فریم‌های خوانش
print("\nترجمه فریم‌های خوانش (با توقف در کدون پایان):")
for frame in range(3):
    protein = genome_segment_record.seq[frame:].translate(table=1, to_stop=True)
    print(f"  فریم {frame}: {protein}")

# 3. شناسایی ORFهای بلندتر از 100 نوکلئوتید
print("\nشناسایی ORFهای بالقوه (طول پروتئین حداقل 10 آمینواسید):")
min_protein_len = 10
standard_table = CodonTable.unambiguous_dna_by_id[1] # جدول استاندارد

orfs_found = []
for strand_idx, dna_strand in enumerate([genome_segment_record.seq, genome_segment_record.seq.reverse_complement()]):
    strand_name = "Forward" if strand_idx == 0 else "Reverse"
    for frame_idx in range(3):
        translated_protein = dna_strand[frame_idx:].translate(table=standard_table, to_stop=False)
        current_orf_start = -1
        for i, aa in enumerate(translated_protein):
            if aa == 'M' and current_orf_start == -1:
                current_orf_start = i * 3 + frame_idx
            elif aa == '*' and current_orf_start != -1:
                orf_protein_seq = translated_protein[current_orf_start // 3 - frame_idx // 3 : i]
                if len(orf_protein_seq) >= min_protein_len:
                    orfs_found.append({
                        "protein_seq": orf_protein_seq,
                        "dna_start": current_orf_start,
                        "dna_end": i * 3 + frame_idx + 3, # شامل کدون پایان
                        "strand": strand_name,
                        "frame": frame_idx
                    })
                current_orf_start = -1
        # اگر ORF به انتهای توالی رسید بدون کدون پایان
        if current_orf_start != -1:
            orf_protein_seq = translated_protein[current_orf_start // 3 - frame_idx // 3 :]
            if len(orf_protein_seq) >= min_protein_len:
                orfs_found.append({
                    "protein_seq": orf_protein_seq,
                    "dna_start": current_orf_start,
                    "dna_end": len(dna_strand),
                    "strand": strand_name,
                    "frame": frame_idx
                })

for orf in orfs_found:
    print(f"  ORF یافت شد: پروتئین '{orf['protein_seq']}', طول: {len(orf['protein_seq'])}, "
          f"موقعیت DNA: {orf['dna_start']}-{orf['dna_end']}, جهت: {orf['strand']}, فریم: {orf['frame']}")

این مثال‌ها نشان می‌دهند که چگونه با ترکیب مفاهیم `Biopython.Seq`، `SeqRecord` و دیگر ابزارهای Biopython می‌توان به تحلیل‌های پیچیده بیولوژیکی پرداخت و insights ارزشمندی را از داده‌های توالی استخراج کرد.

نتیجه‌گیری

در این مقاله جامع، به بررسی عمیق ماژول `Biopython.Seq` و ابزارهای مرتبط آن پرداختیم که نقش محوری در مدیریت و تحلیل توالی‌های DNA و RNA در بیوانفورماتیک نوین ایفا می‌کنند. از مبانی ساخت و دستکاری توالی‌ها با اشیاء `Seq` و `MutableSeq` گرفته تا عملیات‌های بیولوژیکی کلیدی مانند مکمل‌سازی، ترانسکریپشن و ترجمه، و در نهایت نحوه مدیریت ابرداده‌ها با `SeqRecord` و `SeqFeature`، تمام جوانب اصلی مورد بررسی قرار گرفتند.

همچنین، با چگونگی خواندن و نوشتن توالی‌ها به فرمت‌های رایج FASTA و GenBank با استفاده از `Bio.SeqIO` آشنا شدیم و تفاوت‌های اساسی این فرمت‌ها را درک کردیم. در بخش سناریوهای پیشرفته، نمونه‌هایی از کاربردهای عملی مانند شناسایی تفاوت‌های توالی (SNP ها و ایندل‌ها)، ادغام Biopython در پایپ‌لاین‌های بیوانفورماتیکی و تعامل با ابزارهای خارجی مانند BLAST را مشاهده کردیم.

قدرت `Biopython.Seq` تنها در توانایی آن برای انجام عملیات بر روی توالی‌ها نیست، بلکه در ارائه یک چارچوب منسجم و شیءگرا برای نمایش و دستکاری داده‌های بیولوژیکی است که امکان توسعه تحلیل‌های پیچیده و خودکار را فراهم می‌کند. با درک و به‌کارگیری مفاهیم و کدهای ارائه شده در این مقاله، شما اکنون مجهز به دانش و ابزارهایی هستید که می‌توانید پروژه‌های بیوانفورماتیکی خود را با کارایی و دقت بیشتری پیش ببرید.

Biopython یک کتابخانه در حال تکامل است و همواره با قابلیت‌های جدیدی به‌روز می‌شود. توصیه می‌شود برای کسب اطلاعات بیشتر و کاوش در امکانات وسیع‌تر آن، به مستندات رسمی Biopython مراجعه کرده و با تمرین و تجربه، مهارت‌های خود را در این زمینه توسعه دهید. مدیریت توالی‌های DNA و RNA با Biopython.Seq تنها آغاز راهی پربار در دنیای بیوانفورماتیک است که بی‌پایان پتانسیل برای اکتشاف و نوآوری دارد.