Khóa luận Phát hiện marker microsatellite từ cơ sở dữ liệu trình tự EST (Expressed Sequence Tags) của cây xoài (Mangifera indica)

0 trình tự) trong một lần chạy chƣơng trình và không có sự thống kê kết quả tìm kiếm 3.3.3.2. Công cụ MISA (Microsatellite Identification Tool) Công cụ này cho phép sự xác định và định vị microsatellite cũng nhƣ các microsatellite ghép (compound microsatellite) mà bị ngắt quãng bởi một số base nhất định. Công cụ này đƣợc viết từ ngôn ngữ lập trình Perl và có thể tải về máy sử dụng từ địa chỉ internet sau gatersleben.de/misa/misa.html Cú pháp misa.pl với là file chứa trình tự ở định dạng FASTA Công cụ này cần thêm một file chứa các thông số cho quá trình xác định microsatellite, file này có định dạng “misa.ini”. 46 File misa.ini có cấu trúc nhƣ sau Hình 3.8. File misa.ini Trong đó hàng đầu tiên xác định dạng và số lần lặp lại tối thiểu của microsatellite. Hàng thứ hai là số nucleotide tối đa chèn vào giữa hai microsatellite. Kết quả của việc tìm kiếm microsatellite sẽ đƣợc lƣu trong hai file File “.misa” lƣu các giá trị sự định vị và sự xác định các microsatellite File “.statistics” thống kê kết quả của quá trình tìm kiếm 3.3.4. Xác định vùng bảo tồn Vùng bảo tồn là vùng ở hai bên trình tự microsatellite, những trình tự này thƣờng giống nhau đối với mỗi một dạng microsatellite. Vùng bảo tồn rất quan trọng trong phân tích microsatellite bởi đây chính là cơ sở cho việc thiết kế primer. Để xác định đƣợc vùng bảo tồn, chúng tôi tiến hành việc sắp gióng cột (alignment) các trình tự. Sắp gióng cột dựa vào microsatellite, chúng tôi lấy trình tự microsatellite làm điểm tập trung, sau đó tạo trình tự bảo tồn (consensus sequence) đối với mỗi dạng microsatellite bằng công cụ “CAP Contig Assembly Program” có trong phần mềm BioEdit. 47 Hình 3.9. Sắp gióng cột trình tự 3.3.5. Thiết kế primer Primer là những đoạn nucleotide ngắn, bắt cặp bổ sung với đầu 5‟ hay đầu 3‟ của mạch DNA khuôn mẫu. Primer đƣợc thiết kế dựa vào vùng trình tự đã đƣợc biết, nằm ở hai đầu của đoạn gen cần khuếch đại. Thông số quyết định sự thành công của phản ứng PCR là việc thiết kế primer. Một primer đƣợc thiết kế không tốt có thể sẽ cho kết quả ít hay không có sản phẩm do sự khuếch đại không chuyên biệt và/hay sự hình thành cấu trúc thứ cấp, các cấu trúc này sẽ cạnh tranh và ngăn chặn sự tạo thành sản phẩm mong muốn. Việc thiết kế và chọn lựa primer phải thỏa mãn một số yêu cầu sau - Chiều dài primer: chiều dài tốt là 18 đến 24 base. Chiều dài này đủ dài để đảm ảo tính chuyên biệt và đủ ngắn để primer bám vào mạch mẫu dễ dàng ở nhiệt độ bắt cặp. - Nhiệt độ nóng chảy (Tm) là nhiệt độ mà một nửa sợi đôi DNA tách ra trở thành sợi đơn và cho biết tính ổn định của sợi đôi. Thành phần (G+C) trong DNA cao sẽ dẫn tới nhiệt độ Tm cao vì liên kết H trong DNA cao hơn. Có nhiều công thức tính T, hai trong những công thức đƣợc nhiều ngƣời sử dụng là Tm = 59.9 + 0.41 * (%GC) – 675 / chiều dài 48 Tm = 2 (A+C) + 4 (G+C) (công thức Wallace) Primer với nhiệt độ nóng chảy trong phạm vi 55 oC đến 72 oC thƣờng cho kết quả tốt nhất. - Tính chuyên biệt: primer phải đƣợc lựa chọn sao cho chỉ có một trình tự duy nhất trong DNA mẫu đƣợc khuếch đại. Vì Taq polymerase có hoạt tính trong một phạm vi nhiệt độ rộng, sự kéo dài primer sẽ xảy ra ở nhiệt độ thấp hơn nhiệt độ bắt cặp. Nếu sự thay đổi nhiệt độ quá chậm sự không chuyên biệt sẽ xảy ra và enzyme sẽ xúc tác sự kéo dài nếu có một sự tƣơng đồng ngắn ở đầu 3‟. - Thành phần base: ảnh hƣởng đến độ đặc hiệu của quá trình bắt cặp, nhiệt độ nóng chảy, nhiệt độ bắt cặp và sự ổn định của cấu trúc phân tử. Các base đƣợc sắp xếp ngẫu nhiên thì thích hợp hơn là những vùng (A+T) dài hay là những vùng giàu (G+C). Thành phần (G+C) trung bình khoảng từ 50% đến 60% sẽ cho nhiệt độ nóng chảy, nhiệt độ bắt cặp thích hợp trong một phản ứng PCR bình thƣờng. - Trình tự primer đầu 3‟ quyết định tính chuyên biệt và tính tƣơng thích của phản ứng PCR. Các trình tự đầu 3‟ không nên có + G hay C không có nhiều hơn 3 tại vị trí này vì sẽ làm cho primer bắt cặp không chuyên biệt. + 3‟ thymidine, nó làm cho sự bắt cặp sai dễ xảy ra hơn so với những nucleotide khác. Các cặp primer nên đƣợc kiểm tra sự bổ sung ở đầu 3‟ vì nó thƣờng dẫn đến sự hình thành cấu trúc thứ cấp. - Trình tự primer đầu 5‟: các base ở đầu 5‟ ít quyết định đến sự bắt cặp của primer. Do đó có thể thêm các yếu tố trình tự nhƣ vị trí giới hạn (restriction site) - những vị trí này có thể giống nhau hoặc có cùng đầu dính (đầu bằng) với enzyme giới hạn trong MCS (Multiple Cloning Site) của vector chọn để dòng hóa gen quan tâm. NcoI (CCATGG) hay NdeI (CATATG) thƣờng đƣợc sử dụng vì có thể tạo ra codon mở đầu ATG. - Cấu trúc thứ cấp: nếu sự bắt cặp giữa forward primer với reverse primer (hình thành dimer, hetero-dimer), forward primer với forward 49 primer hay reverse primer với reverse primer (self-dimer, homo- dimer), hay primer tự tạo thành cấu trúc hairpin (kẹp tóc) xảy ra nhiều hơn so với sự bắt cặp của primer với DNA mẫu thì hiệu quả nhân bản của phản ứng PCR sẽ giảm một cách rõ rệt. Nhƣ vậy, nên tránh những trƣờng hợp này. 3.3.5.1. Primer3 Chƣơng trình Primer3 là một trong những chƣơng trình thiết kế primer miễn phí đƣợc biết đến nhiều nhất. Chƣơng trình này đƣợc tạo ra bởi các nhà khoa học thuộc Viện Nghiên cứu Sinh Y học Whitehead và Trung Tâm Nghiên cứu Genome của MIT (Whitehead Institute for Biomedical Research and MIT Center for Genome Research). Chƣơng trình này có thể đƣợc sử dụng hoàn toàn miễn phí tại địa chỉ internet bin/primer/primer3_www.cgi Chƣơng trình Primer3 thiết kế primer cho một trình tự DNA đƣa vào, thỏa mãn nhiều tùy chọn khác nhau, chủ yếu là các điều kiện về %GC, về nhiệt độ bắt cặp của primer, kích thƣớc sản phẩm… Đây là chƣơng trình lớn với gần 100 tùy chọn khác nhau tƣơng ứng với các điều kiện mà primer đƣợc tạo ra phải thỏa mãn. Có một thuận lợi là hầu hết các tùy chọn này đều có giá trị mặc định của nó, và ngƣời dùng có thể không thay đổi các thông số này nếu họ không có nhu cầu đặc biệt. Khi muốn thiết kế primer cho một gen hay một đoạn trình tự nào đó, trƣớc tiên ngƣời sử dụng phải đƣa đoạn trình tự DNA của mình vào chƣơng trình và xác định các thông số về trình tự primer, nhiệt độ nóng chảy cua primer, kích thƣớc sản phẩm cần… Các thông số chi tiết khác có thể giữ mặc định. 50 Hình 3.10. Chƣơng trình Primer3 Dựa trên những yêu cầu cho việc thiết kế primer đã xác định nhƣ trên chúng tôi thiết lập các thông số cho chƣơng trình Primer3 nhƣ sau + Targets(mục tiêu): m, n với m là vị trí bắt đầu có microsatllite và n là chiều dài microsatellite +Primer size (kích thƣớc primer) Min: 18; Max: 24 + Primer Tm (nhiệt độ nóng chảy của pimer) Min: 55.0; Max: 72.0; Max Tm difference: 2.0 + Primer %GC Min: 50.0; Max: 60% + Các thông số khác vẫn giữ mặc định. 3.3.5.2. Chƣơng trình Perl ssrfinder_1_0 Đây là một chƣơng trình của tác giả Steven Schroeder thuộc trƣờng Đại học Missouri – Michigan. Chƣơng trình gồm 6 Perl scripts có chức năng xác định SSR và thiết kế primer thích hợp cho mỗi SSR tìm đƣợc 51 - 1_ssr_repeat_finder.pl: tìm SSR, lấy ra trình tự SSR và vùng flanking cho những phân tích sau. - 2_ssr_primer_designer.pl: thiết kế primer mà mục tiêu là khuếch đại vùng trình tự chứa SSR. - 3_ssr_primer_rep_check.pl: sàng lọc lại các primer đã thiết kế để loại bỏ những primer có chứa trình tự lặp lại - 4_ssr_primer_blast.pl: so sánh các primer đã thiết kế với cơ sở dữ liệu primer - 5_ssr_order_filter.pl: tạo 1 file chỉ chứa SSR mà có primer duy nhất - 6_ssr_primer_formatter.pl: tạo 1 file chỉ chứa SSR có primer duy nhất – file này đƣợc tạo đơn giản chỉ chứa những thông tin cần thiết cho việc chọn lựa primer. Yêu cầu: vì chƣơng trình này đƣợc viết cho hệ điều hành Unix hay Linux nên cần phải thực hiện sửa đổi một số lệnh lập trình cơ bản để có thể chạy trên môi trƣờng Window. Chƣơng trình cần sự kết hợp với 3 phần mềm khác là Primer3, blastall và formatdb để thực thi. Ba phần mềm này có thể tải hoàn toàn miễn phí (có phiên bản dành cho Window) từ trang Primer3 và trang Blast của NCBI Các thông số của chƣơng trình thiết kế primer đều đƣợc mặc định nhƣ sau TARGET= m, n với m là vị trí bắt đầu có microsatllite và n là chiều dài microsatellite (mục tiêu) PRIMER_PRODUCT_SIZE_RANGE=80-160 80-240 80-300 (kích thƣớc sản phẩm) PRIMER_OPT_SIZE=24 (kích thƣớc tối ƣu của primer) PRIMER_MIN_SIZE=20 (kích thƣớc tối thiểu của primer) PRIMER_MAX_SIZE=28 (kích thƣớc tối đa của primer) PRIMER_OPT_TM=63 (nhiệt độ nóng chảy tối ƣu của primer) PRIMER_MIN_TM=60 (nhiệt độ nóng chảy tối thiểu của primer) PRIMER_MAX_TM=65 (nhiệt độ nóng chảy tối đa của primer) 52 PRIMER_MAX_DIFF_TM=1 (độ chênh lệch nhiệt độ nóng chảy tối đa) Với các thông số mặc định trên chƣơng trình hoàn toàn có thể sử dụng cho mục tiêu của đề tài. 53 Phần 4 KẾT QUẢ VÀ THẢO LUẬN 4.1. Thu nhận trình tự EST của cây xoài Khi sử dụng từ khóa và phƣơng pháp nhƣ mục 3.3.1.2 chúng tôi đã tải đƣợc toàn bộ 15966 trình tự EST của cây xoài. Các trình tự này ở định dạng FASTA, mỗi trình tự có cấu trúc nhƣ hình 4.1 Hình 4.1 Trình tự EST ở định dạng FASTA Đây là kiểu định dạng phổ biến trong cơ sở dữ liệu của GenBank gồm 2 phần chủ yếu. Phần 1 bắt đầu bằng một dấu “>”, theo sau là các thông tin về trình tự đó nhƣ tên trình tự, số gi, accession number…Phần 2 là các nucleotide của trình tự. Định dạng FASTA có ƣu điểm là chỉ chứa trình tự và những thông tin thiết yếu về trình tự đó. Hơn nữa trong một file trình tự ở định dạng FASTA có thể có nhiều trình tự so với dạng Plain format chỉ chứa một trình tự. Mặc khác việc lựa chọn kiểu định dạng trình tự để tải về ở dạng FASTA là sự thuận lợi cho các nghiên cứu sau vì phần lớn các chƣơng trình, phần mềm của tin sinh học đều sử dụng định dạng FASTA. 54 Do mục tiêu của đề tài là tìm kiếm tất cả các microsatellite có trong toàn bộ nguồn dữ liệu EST của cây xoài nên việc tải tất cả 15966 trình tự thuộc cùng một file là điều hợp lý, dễ dàng cho các công việc phân tích sau này. 4.2. Sắp xếp các trình tự Sau khi chạy chƣơng trình est_trimmer.pl với các thông số nhƣ mục 3.3.2 đã trình bày, kết quả đã đƣợc lƣu trong 2 file mango.results: chứa tất cả các trình tự thỏa mãn các thông số của bƣớc kiểm tra với est_trimmer.pl mango.txt: file này ghi nhận tất cả các tiến trình thực thi trên mỗi trình tự không đạt yêu cầu đề ra nhƣ loại bỏ trình tự có kích thƣớc bé hơn 100 bp, loại bỏ đuôi poly A, poly T… Hình 4.2. Tiến trình thực thi của est_trimmer.pl 4.3. Kết quả tìm kiếm microsatellite 4.3.1. Công cụ SSRIT SSRIT cho phép tìm kiếm một lần là 500 trình tự và tốn khoảng 5 phút. Chính vì công cụ SSRIT không có khả năng tìm kiếm SSR với một số lƣợng quá lớn trình tự 55 (15966) nên phải chạy chƣơng trình nhiều lần, tốn nhiều thời gian và phải lặp lại tiến trình. Hơn nữa bảng kết quả của SSRIT không lƣu thành file để truy cập, phải xem kết quả trực tuyến, không có sự thống kê về các dạng microsatellite. Do đó chúng tôi xác định rằng công cụ này không phù hợp với mục đích nghiên cứu của đề tài. 4.3.2. Công cụ MISA Thực thi chƣơng trình MISA cho kết quả rất nhanh và có 2 file đƣợc tạo thành mango.fasta.misa: chứa các thông tin về tên trình tự, dạng SSR, kích thƣớc, vị trí bắt đầu và kết thúc của SSR… Hình 4.3 Nội dung file mango.fasta.misa mango.fasta.stastistics: là file thống kê kết quả tìm kiếm microsatellite Dựa vào file mango.fasta.statistics này chúng tôi xác định đƣợc số lƣợng trình tự không đáp ứng yêu cầu mà khi chạy với est_trimmer.pl đã chƣa thống kê đƣợc: ban đầu chúng tôi có 15966 trình tự, sau khi chạy est_trimmer thì có 231 trình tự không thỏa mãn yêu cầu do đó tổng số trình tự đƣợc kiểm tra là 15735. Kết quả của việc tìm kiếm microsatellite đƣợc liệt kê ở bảng 4.1, số lƣợng SSR đƣợc xác định chiếm tỷ lệ 1.7% tổng số trình tự EST đƣợc kiểm tra, phân bố tƣơng đối nhiều trong bộ gen. Với nghiên cứu trên đối tƣợng là cây nho (Scott và ctv, 2000), 56 SSR chiếm 2.5% và trên cây bông vải (Qureshi, 2004) SSR chiếm 1.34% tổng số trình tự EST. Bảng 4.1. Kết quả tìm kiếm microsatellite Kết quả Số lƣợng Tổng số trình tự kiểm tra 15735 Tổng kích thƣớc của trình tự đã kiểm tra (bp) 7400551 Tổng số SSR đƣợc xác định 267 Tổng số trình tự có SSR 265 Số trình tự có nhiều hơn 1 SSR 2 Số SSR hiện diện ở dạng SSR ghép 1 Sự phân bố các dạng lặp lại có số lƣợng khác nhau theo bảng 4.2. Bảng 4.2. Sự phân bố các dạng lặp lại của SSR Đơn vị lặp lại Số SSR Tỷ lệ % 2 (dinucleotide) 11 4.12 3 (trinucleotide) 255 95.51 4 (tetranucleotide) 1 0.37 Bảng trên cho thấy phần lớn các microsatellite của cây xoài là dạng trinucleotide (95.51%). Các dạng khác xuất hiện rất ít, điều này tƣơng tự đối với microsatllite ghép (compound microsatellite). Trong các nghiên cứu khác trên đối tƣợng cây nho dạng trinucleotide là chỷ yếu chiếm tỷ lệ là 62.90%. Trên cây mía dạng dinucleotide có tỷ lệ 38.19%, trong khi dạng trinucleotide có tỷ lệ 31.49%. Sự phân bố các dạng lặp lại của microsatellite có sự khác biệt giữa các đối tƣợng nghiên cứu khác nhau. Ƣu điểm của MISA so với các phƣơng pháp khác là kết quả tìm kiếm của MISA đƣợc lƣu thành file, có sự thống kê các dạng microsatellite, chƣơng trình thực thi nhanh, máy tính không cần nối mạng internet chỉ cần tải mã (code) của chƣơng trình về máy tính sau đó có thể áp dụng đơn giản và hiệu quả. 57 Nhƣợc điểm của chƣơng trình này là không có chức năng phân loại các trình tự theo dạng microsatellite đã tìm kiếm. Do đó chúng tôi phải tiến hành phân nhóm trình tự bằng cách thủ công dựa theo các dạng microsatellite đã xác định đƣợc Kết quả đạt đƣợc cho thấy Có tất cả 31 loại SSR trong đó dạng dinucleotide có 5 loại, dạng trinucleotide có 25 loại và tetranucleotide có 1 loại. Bảng 4.3. Các loại SSR Dạng dinucleotide Dạng trinucleotide Dạng tetranucleotide AT; CA; CT; TC; GA AAC; AAG; ACC; AGA; AGC; AGT; ATC; ATT; CAA; CAG; CAT; CCA; CGA; CTT; GAA; GCC; GCT; GGA; GGT; GTG; TCA; TCT; TGA; TGG; TTC TAAT Hình 4.4 Các file trình tự sau khi phân nhóm 58 Vì các dạng dinucleotide và tetranucleotide chiếm tỷ lệ và tần suất xuất hiện thấp nên chúng tôi chỉ tiếp tục nghiên cứu đối với dạng trinucleotide. 4.4. Xác định vùng bảo tồn Để kết quả align (sắp gióng cột) và xác định vùng bảo tồn đƣợc chính xác, chúng tôi loại bỏ các loại microsatellite có tần suất xuất hiện thấp (<2.5%) Sau quá trình sàng lọc chúng tôi còn lại 7 loại microsatellite Bảng 4.4. Các loại microsatellite nghiên cứu Số thự tự Loại microsatellite Số trình tự Tần suất xuất hiện 1 AGA 31 11.61% 2 CAA 14 5.24% 3 CAT 24 8.99% 4 CCA 73 27.34% 5 TCA 7 2.62% 6 TCT 32 11.99% 7 TGA 16 5.99% Thực hiện việc sắp gióng cột và tạo trình tự bảo tồn dựa trên công cụ CAP contig assembly program chúng tôi có đƣợc 6 loại microsatellite, microsatellite AGA sau khi thực hiện 2 tiến trình nói trên cho thấy sự khác biệt về trình tự ở vùng flanking nên không thể thiết kế đƣợc primer chung cho tất cả 31 trình tự. 59 Hình 4.5. Xác định vùng bảo tồn của microsatellite CAA 4.5. Thiết kế primer đối với 6 microsatellite 4.5.1. Primer3 Việc thiết kế primer và chọn lựa các thông số phải thỏa các yêu cầu đã nêu trong mục 3.3.5.1. Kết quả của chƣơng trình Primer3 của 6 microsatellite Hình 4.6. Kết quả thiết kế primer của microsatellite TCA Chúng tôi đã thiết kế đƣợc 6 cặp primer cho 6 loại microsatellite. Theo bảng 4.5, chúng tôi liệt kê loại SSR, forward primer, reverve primer và nhiệt độ nóng chảy (Tm) Bảng 4.5 Kết quả thiết kế primer từ chƣơng trình Primer3 Loại SSR Forward primer (5’ - 3’) Reverse primer (5’ - 3’) CAA GTCTCAAACCCCAAATCC TGTCCACCAGAACCAGAG CAT ACCGCCAGTCTATGGAACAC TACTTAGCGGCCATCAAACC CCA TTCCAAGAAGCACTGTGTGG TTGTGGTTGAGGTTGCTCAG TCA GCATATGCAACTCCACCAG GTTGTTGTGGTGGTCCAA TCT TGCTGCAATCGAACCAAG GCAGGAGGCACAAATAGTAGTC 60 TGA AGCAGTGAATCCGGTGATCT GGTTTCAGCGCCTTATCAAC 4.5.2. Chƣơng trình Perl script ssrfinder_1_0 Lần lƣợt thực thi từng file từ 1 đến 6 của gói công cụ này, chúng tôi có đƣợc kết quả sau  1_ssr_repeat_finder.pl có 3 file kết quả new_id20060715.txt: liệt kê các số truy cập, tên trình tự có SSR labdout20060715.txt: liệt kê dạng SSR, vị trí và chiều dài SSR, trình tự ssrout20060715.txt: liệt kê SSR, trình tự có chứa SSR và vùng flanking  2_ssr_primer_designer.pl có 3 file kết quả primerin.txt: liệt kê các thông số mặc định cho quy trình thiết kế primer (xem mục 3.3.5.3) raw_primer20060715.out: liệt kê các kết quả tính toán cho mỗi primer đƣợc thiết kế primer_result20060715.txt: liệt kê SSR, trình tự, cặp primer đã thiết kế, nhiệt độ nóng chảy và chiều dài sản phẩm của primer. Hình 4.7 Nội dung file primer_result20060715.txt  3_ssr_primer_check.pl tạo file rescreened20060715.txt có nội dung giống file primer_result20060715.txt vì không có primer đƣợc thiết kế có nhiều trình tự lặp lại nên bƣớc kiểm tra này vẫn giữ nguyên giá trị cũ.  4_ssr_primer_blast20060715.pl cũng cho kết quả giống primer_result20060715.txt vì cơ sở dữ liệu primer hiện tại không có.  5_ssr_order_filter.pl: chƣơng trình này tạo file filter20060715.txt có chứa dạng SSR, trình tự primer, và các thông số của primer khác  6_ssr_order_formatter tạo file cuối cùng là order20060715.txt, file này chỉ chứa các thông tin cần thiết cho việc thiết kế mồi nhƣ dạng microsatellite, forward primer và reverse primer. 61 * Kết quả của chƣơng trình thiết kế primer Hình 4.8. Kết quả thiết kế primer Trong hình 4.8 cột A là loại SSR, cột B là forward primer, cột C và E là nhiệt độ nóng chảy, cột D là reverse primer, cột F là chiều dài sản phẩm tạo thành. Từ kết quả thiết kế primer thu đƣợc từ 2 chƣơng trình khác nhau, chúng tôi nhận thấy các primer trên đều thích hợp cho phản ứng PCR tùy vào yêu cầu, mục đích thực nghiệm. Primer ở 2 chƣơng trình có kết quả khác nhau là do việc lựa chọn các thông số cho tiến trình thiết kế mồi khác nhau. 62 Phần 5 KẾT LUẬN VÀ ĐỀ NGHỊ 5.1. Kết luận Từ quá trình thực hiện đề tài chúng tôi rút ra kết luận nhƣ sau Xây dựng phƣơng pháp phát hiện microsatellite đối với cây xoài từ nguồn cơ sở dữ liệu EST hiện có, đó là sử dụng các chƣơng trình và công cụ est_trimmer.pl, misa.pl, CAP contig assembly program, Primer3 hay gói công cụ ssrfinder_1_0. Cơ sở dữ liệu EST của NCBI Sử dụng công cụ tải dữ liệu Truy cập của NCBI Thu nhận Xử lý với chƣơng trình est_trimmer.pl Sắp xếp các trình tự EST Xử lý với chƣơng trình misa.pl Xác định microsatellite Xử lý với phần mềm BioEdit Align, xác định vùng bảo tồn Xử lý với chƣơng trình Primer3, Thiết kế primer gói công cụ ssrfinder_1_0 Hình 5.1. Sơ đồ phƣơng pháp thực hiện Trình tự EST của cây xoài Các trình tự EST đã đƣợc sắp xếp Các trình tự EST có chứa microsatellite Align và xác định vùng bảo tồn Primer đã thiết kế 63 Kết quả đạt đƣợc Đã tải đƣợc 15966 trình tự EST của cây xoài Xác định và định vị đƣợc 267 microsatellite trong đó dạng trinucleotide chiếm tỉ lệ cao nhất (95.51%) Xác định đƣợc vùng bảo tồn của 6 loại microsatellite và thiết kế đƣợc primer cho 6 loại microsatellite bao gồm các loại sau: CAA, CAT, CCA, TCA, TCT và TGA. 5.2 Đề nghị Đề tài nên đƣợc tiếp tục phân tích bằng thực nghiệm để kiểm tra hiệu quả của việc thiết kế primer. Viết và sử dụng chƣơng trình hay phần mềm có khả năng phân loại các trình tự đã xác định SSR. Áp dụng quy trình xác định microsatellite cho các đối tƣợng thí nghiệm khác nếu đối tƣợng thí nghiệm đó có sẵn nguồn dữ liệu từ các cơ sở dữ liệu. Xây dựng cơ sở dữ liệu trình tự, microsatellite và các primer đã thiết kế. Tăng cƣờng việc đào tạo, phát triển ngành tin sinh học vì đó là một ngành hỗ trợ đắc lực cho công tác nghiên cứu. 64 Phần 6 TÀI LIỆU THAM KHẢO TIẾNG VIỆT 1. Phạm Văn Duệ, 2005. Giáo trình kỹ thuật trồng cây ăn quả. Nhà xuất bản Hà Nội. Trang 105. 2. Phạm Thị Hƣơng, Trần Thế Tục, Nguyễn Quang Thạch, 2003. Cây xoài và những điều cần biết. Nhà xuất bản Nông Nghiệp Hà Nội. 95 trang. 3. Trần Thế Tục, 1994. Kỹ thuật trồng và chăm sóc Xoài Na Hồng Xiêm. Nhà xuất bản Nông Nghiệp Hà Nội. Trang 25. 4. Lê Minh Trung, Quốc Bình, 2002. Ngôn ngữ lập trình Perl cho người mới học. Nhà xuất bản Thống kê. 426 trang. 5. Trần Văn Lăng, 2003. Một số tổng quan về sinh tin học. Phân viện Công nghệ Thông tin tại TP Hồ Chí Minh. 6. Trần Linh Thƣớc, 2004. Thực tập bioinformatics. Đại học Khoa học Tự nhiên TP Hồ Chí Minh. 7. Lƣu Phúc Lợi, 2005. Giáo trình sinh tin học ứng dụng. Bộ môn Công nghệ Sinh học – Đại học Nông Lâm TP Hồ Chí Minh. TIẾNG NƢỚC NGOÀI 8. Andreas D. Baxevanis, B. F. Francis Ouellette, 2005. Bioinformatics a practical guide to the analysis of genes and proteins. 3 rd edition, A John Wiley & Sons, Inc., Publication. 9. James Tisdall, 2001. Beginning Perl for Bioinformatics. 1st edition, O‟Reilly, Publication. 10. Lorenzo Cerutti, 2002. EST clustering. Swiss Institute of Bioinformatics 11. Scott K. D., 2001. Microsatellites derived from ESTs, and their comparison with those derived by other methods. Centre for plant conservation genetics, Southern Cross University, Lismore, Australia. 12. Scott K. D., P. Eggler, G. Seaton, M. Rossetto, E. M. Ablett, L. S. Lee, R. J. Henry, 2000. Analysis of SSRs derived from grape ESTs. 65 13. Samina N. Qureshi, Sukumar Saha, Ramesh V. Kantety, and J. N. Jenkins, 2004. EST-SSR: A New Class of Genetic Markers in Cotton. The journal of cotton science 8:112 – 123. The cotton foundation. TÀI LIỆU TỪ CÁC TRANG WEB 14. Julian P. Robinson, Stephen A. Harris, 1999. Amplified Fragment Length Polymorphisms and Microsatellites: A phylogenetic perspective. Department of Plant Sciences, University of Oxford, South Parks Road, Oxford OX1 3RB, Great Britain 15. Vincent R. Prezioso. General notes on Primer Design in PCR. Biosystems Laboratory, Brinkmann Instruments Inc., Westbury, New York. 16. 17. 18. 19. 20. 21. 22. 23. 66 PHỤ LỤC 1. Mã chƣơng trình est_trimmer.pl #!/usr/bin/perl –w # Check for arguments. If none display syntax # if (@ARGV == 0) {open (IN,"<$0"); while () {if (/^\#\# (.*)/) {$message .= "$1\n"}}; close (IN); die $message;} # Check if help is required # if ($ARGV[0] =~ /-help/i) {open (IN,"<$0"); while () {if (/^\#\#\#(.*)/) {$message .= "$1\n"}}; close (IN); die $message;} # Open FASTA file # open (IN,"<$ARGV[0]") || die ("\nError: File doesn't exist !\n\n"); # Checking arguments # $arg = @ARGV; $arg > 1 || die ("\nError: No arguments determined !\n\n"); for ($i = 1; $i < $arg; $i++) {if (($amb_n,$amb_win) = ($ARGV[$i] =~ /-amb=(\d+),(\d+)/i)) {$message .= "$i. Check for ambiguous bases (search for $amb_n ambiguous bases in a $amb_win bp window).\n";} elsif (($tr3_b,$tr3_n,$tr3_win) = ($ARGV[$i] =~ /- tr3=([ACGT]),(\d+),(\d+)/i)) {$message .= "$i. Trim 3' end: Remove \"$tr3_b\" tails."; if ($tr3_win > 0) {$message .= " Check for $tr3_n x $tr3_b in a $tr3_win bp window.\n"} else {$message .= "\n"};} elsif (($tr5_b,$tr5_n,$tr5_win) = ($ARGV[$i] =~ /- tr5=([ACGT]),(\d+),(\d+)/i)) {$message .= "$i. Trim 5' end: Remove \"$tr5_b\" tails."; if ($tr5_win > 0) {$message .= " Check for $tr5_n x $tr5_b in a $tr5_win bp window.\n"} else {$message .= "\n"};} elsif (($cut_min,$cut_max) = ($ARGV[$i] =~ /-cut=(\d+),(\d+)/i)) {$message .= "$i. Size restrictions: Size cutoff is $cut_min bp. Restrict sequence size to $cut_max bp.\n";} elsif (($file) = ($ARGV[$i] =~ /-id=(.*)/i)) {} 67 else {die ("\nError: Argument nr. ",++$i," is invalid !\n\n")}}; # Open results & log file# if ($file) {open (OUT,">$file.results") || die ("\nError: Output file name not valid !\n\n");open (LOG,">$file.log")} else {open (OUT,">$ARGV[0].results");open (LOG,">$ARGV[0].log")}; print "\nEST TRIMMER\n===========\n\nPerforming following steps:\n\n"; print LOG "\nLOG FILE OF ALL PERFORMED MODIFICATIONS\n"; print LOG "=======================================\n\n"; print LOG "\nPerforming following steps:\n--------------------------- \n\n"; print $message; print LOG $message; print LOG "\n\nProcessing steps for each sequence:\n-------------------- ---------------\n"; # core # $/ = ">"; while () {next unless (($seqname,$seq) = /(.*?)\n(.*)/s); $seq =~ s/[\d\s>]//g; $seq =~ s/[^ACGTN]/N/gi; # only "acgtn" characters allowed $message = ''; $discard = 'nö'; for ($i = 1; $i < $arg; $i++) {if (($amb_n,$amb_win) = ($ARGV[$i] =~ /-amb=(\d+),(\d+)/i)) {&Trim_5_n; &Trim_3_n;} elsif (($tr3_b,$tr3_n,$tr3_win) = ($ARGV[$i] =~ /- tr3=([ACGT]),(\d+),(\d+)/i)) { &Trim_3_oligoN } elsif (($tr5_b,$tr5_n,$tr5_win) = ($ARGV[$i] =~ /- tr5=([ACGT]),(\d+),(\d+)/i)) {&Trim_5_oligoN } elsif (($cut_min,$cut_max) = ($ARGV[$i] =~ /-cut=(\d+),(\d+)/i)) {if ($cut_max < length $seq) {$seq = substr $seq,0,$cut_max; $message .= "Restrict sequence size to $cut_max bp.\n"}; if ($cut_min and (length $seq < $cut_min)) {$message .= "Discard sequence (below $cut_min bp cutoff).\n"; $discard = "yo;} if ($discard eq 'nö') {$seq = join ("\n",grep($_,split(/(.{70})/,$seq))); print OUT ">$seqname\n$seq\n"}; 68 if ($message ne '') {print LOG "\n>$seqname\n$message"} }; print "\nDONE\n"; close (IN); close (OUT); close (LOG); # subroutines # sub Trim_5_n # trim 5' end until specified window contains less than n ambiguous bases {$check = '0'; while ((length $seq > $amb_win) and !($check eq '1')) {$window = substr $seq,0,$amb_win; if ($window =~ /^(([^N]*N){$amb_n})/i) {$seq =~ s/^($1N*)//i; $message .= "Ambiguous sequence at 5' side: $1.\n" } else {$check = '1'};}; sub Trim_3_n # trim 3' end until specified window contains less than n ambiguous bases {$check = '0'; while ((length $seq > $amb_win) and !($check eq '1')) {$window = substr $seq,-$amb_win; if ($window =~ /.*((N[^N]*){$amb_n})/i) { $seq =~ s/(N*$1)$//i; $message .= "Ambiguous sequence at 3' side: $1.\n" } else {$check = '1'}; }; sub Trim_5_oligoN # remove oligoN stretches from 5' side {$check = '0'; if ($seq =~ s/^($tr5_b+)//i) {$message .= "Remove \"$tr5_b\" tail at 5' side: $1.\n"}; while (!($check eq '1')) { if (length $seq > $tr5_win) {$window = substr $seq,0,$tr5_win} else {$window = $seq}; if ($window =~ /^(.*?($tr5_b){$tr5_n})/i) { $seq =~ s/^($1$tr5_b*)//i; $message .= "Remove \"$tr5_b\" stretch at 5' side: $1.\n" } else {$check = '1'}; }; sub Trim_3_oligoN # remove oligoN stretches from 3'side {$check = '0'; 69 if ($seq =~ s/($tr3_b+)$//i) {$message .= "Remove \"$tr3_b\" tail at 3' side: $1.\n"}; while (!($check eq '1')) {if (length $seq > $tr3_win) {$window = substr $seq,-$tr3_win} else {$window = $seq}; if ($window =~ /.*(($tr3_b){$tr3_n}.*)/i) { $seq =~ s/($tr3_b*$1)$//i; $message .= "Remove \"$tr3_b\" stretch at 3' side: $1.\n"} else {$check = '1'};} }; 2. Mã chƣơng trình misa.pl #!/usr/bin/perl –w #§§§§§ DECLARATION §§§§§# # Check for arguments. If none display syntax # if (@ARGV == 0) {open (IN,"<$0"); while () {if (/^\#\# (.*)/) {$message .= "$1\n"}}; close (IN); die $message;}; # Check if help is required # if ($ARGV[0] =~ /-help/i) {open (IN,"<$0"); while () {if (/^\#\#\#(.*)/) {$message .= "$1\n"}}; close (IN); die $message;}; # Open FASTA file # open (IN,"<$ARGV[0]") || die ("\nError: FASTA file doesn't exist !\n\n"); open (OUT,">$ARGV[0].misa"); print OUT "ID\tSSR nr.\tSSR type\tSSR\tsize\tstart\tend\n"; # Reading arguments # open (SPECS,"misa.ini") || die ("\nError: Specifications file doesn't exist !\n\n"); my %typrep; my $amb = 0; while () { %typrep = $1 =~ /(\d+)/gi if (/^def\S*\s+(.*)/i); if (/^int\S*\s+(\d+)/i) {$amb = $1}}; my @typ = sort { $a $b } keys %typrep; #§§§§§ CORE §§§§§# $/ = ">"; my $max_repeats = 1; #count repeats 70 my $min_repeats = 1000; #count repeats my (%count_motif,%count_class); #count my ($number_sequences,$size_sequences,%ssr_containing_seqs); #stores number and size of all sequences examined my $ssr_in_compound = 0; my ($id,$seq); while () {next unless (($id,$seq) = /(.*?)\n(.*)/s); my ($nr,%start,@order,%end,%motif,%repeats); # store info of all SSRs from each sequence $seq =~ s/[\d\s>]//g; #remove digits, spaces, line breaks,... $id =~ s/^\s*//g; $id =~ s/\s*$//g;$id =~ s/\s/_/g; #replace whitespace with "_" $number_sequences++; $size_sequences += length $seq; for ($i=0; $i < scalar(@typ); $i++) #check each motif class { my $motiflen = $typ[$i]; my $minreps = $typrep{$typ[$i]} - 1; if ($min_repeats > $typrep{$typ[$i]}) {$min_repeats = $typrep{$typ[$i]}}; #count repeats my $search = "(([acgt]{$motiflen})\\2{$minreps,})"; while ( $seq =~ /$search/ig ) #scan whole sequence for that class { my $motif = uc $2; my $redundant; #reject false type motifs [e.g. (TT)6 or (ACAC)5] for ($j = $motiflen - 1; $j > 0; $j--) { my $redmotif = "([ACGT]{$j})\\1{".($motiflen/$j-1)."}"; $redundant = 1 if ( $motif =~ /$redmotif/ )}; next if $redundant; $motif{++$nr} = $motif; my $ssr = uc $1; $repeats{$nr} = length($ssr) / $motiflen; $end{$nr} = pos($seq); $start{$nr} = $end{$nr} - length($ssr) + 1; # count repeats $count_motifs{$motif{$nr}}++; #counts occurrence of individual motifs $motif{$nr}->{$repeats{$nr}}++; #counts occurrence of specific SSR in its appearing repeat $count_class{$typ[$i]}++; #counts occurrence in each motif class if ($max_repeats < $repeats{$nr}) {$max_repeats = $repeats{$nr}};}; }; next if (!$nr); #no SSRs $ssr_containing_seqs{$nr}++; 71 @order = sort { $start{$a} $start{$b} } keys %start; #put SSRs in right order $i = 0; my $count_seq; #counts my ($start,$end,$ssrseq,$ssrtype,$size); while ($i < $nr) { my $space = $amb + 1; if (!$order[$i+1]) #last or only SSR {$count_seq++; my $motiflen = length ($motif{$order[$i]}); $ssrtype = "p".$motiflen; $ssrseq = "($motif{$order[$i]})$repeats{$order[$i]}"; $start = $start{$order[$i]}; $end = $end{$order[$i++]}; next}; if (($start{$order[$i+1]} - $end{$order[$i]}) > $space) { $count_seq++; my $motiflen = length ($motif{$order[$i]}); $ssrtype = "p".$motiflen; $ssrseq = "($motif{$order[$i]})$repeats{$order[$i]}"; $start = $start{$order[$i]}; $end = $end{$order[$i++]}; next }; my ($interssr); if (($start{$order[$i+1]} - $end{$order[$i]}) < 1) { $count_seq++; $ssr_in_compound++; $ssrtype = 'c*'; $ssrseq = "($motif{$order[$i]})$repeats{$order[$i]}($motif{$order[$i+1]})$repeats{$or der[$i+1]}*"; $start = $start{$order[$i]}; $end = $end{$order[$i+1]} } else {$count_seq++; $ssr_in_compound++; $interssr = lc substr($seq,$end{$order[$i]},($start{$order[$i+1]} - $end{$order[$i]}) - 1); $ssrtype = 'c'; $ssrseq = "($motif{$order[$i]})$repeats{$order[$i]}$interssr($motif{$order[$i+1]})$re peats{$order[$i+1]}"; $start = $start{$order[$i]}; $end = $end{$order[$i+1]}; #$space -= length $interssr }; while ($order[++$i + 1] and (($start{$order[$i+1]} - $end{$order[$i]}) <= $space)) { if (($start{$order[$i+1]} - $end{$order[$i]}) < 1) 72 { $ssr_in_compound++; $ssrseq .= "($motif{$order[$i+1]})$repeats{$order[$i+1]}*"; $ssrtype = 'c*'; $end = $end{$order[$i+1]}} else { $ssr_in_compound++; $interssr = lc substr($seq,$end{$order[$i]},($start{$order[$i+1]} - $end{$order[$i]}) - 1); $ssrseq .= "$interssr($motif{$order[$i+1]})$repeats{$order[$i+1]}"; $end = $end{$order[$i+1]}; #$space -= length $interssr}}; $i++;} continue {print OUT "$id\t$count_seq\t$ssrtype\t$ssrseq\t",($end - $start + 1),"\t$start\t$end\n"};}; close (OUT); open (OUT,">$ARGV[0].statistics"); #§§§§§ INFO §§§§§# #§§§ Specifications §§§# print OUT "Specifications\n==============\n\nSequence source file: \"$ARGV[0]\"\n\nDefinement of microsatellites (unit size / minimum number of repeats):\n"; for ($i = 0; $i < scalar (@typ); $i++) {print OUT "($typ[$i]/$typrep{$typ[$i]}) "};print OUT "\n"; if ($amb > 0) {print OUT "\nMaximal number of bases interrupting 2 SSRs in a compound microsatellite: $amb\n"}; print OUT "\n\n\n"; #§§§ OCCURRENCE OF SSRs §§§# #small calculations my @ssr_containing_seqs = values %ssr_containing_seqs; my $ssr_containing_seqs = 0; for ($i = 0; $i < scalar (@ssr_containing_seqs); $i++) {$ssr_containing_seqs += $ssr_containing_seqs[$i]}; my @count_motifs = sort {length ($a) length ($b) || $a cmp $b } keys %count_motifs; my @count_class = sort { $a $b } keys %count_class; for ($i = 0; $i < scalar (@count_class); $i++) {$total += $count_class{$count_class[$i]}}; #§§§ Overview §§§# print OUT "RESULTS OF MICROSATELLITE SEARCH\n================================\n\n"; 73 print OUT "Total number of sequences examined: $number_sequences\n"; print OUT "Total size of examined sequences (bp): $size_sequences\n"; print OUT "Total number of identified SSRs: $total\n"; print OUT "Number of SSR containing sequences: $ssr_containing_seqs\n"; print OUT "Number of sequences containing more than 1 SSR: ",$ssr_containing_seqs - ($ssr_containing_seqs{1} || 0),"\n"; print OUT "Number of SSRs present in compound formation: $ssr_in_compound\n\n\n"; #§§§ Frequency of SSR classes §§§# print OUT "Distribution to different repeat type classes\n----------------- ----------------------------\n\n"; print OUT "Unit size\tNumber of SSRs\n"; my $total = undef; for ($i = 0; $i < scalar (@count_class); $i++) {print OUT "$count_class[$i]\t$count_class{$count_class[$i]}\n"}; print OUT "\n"; #§§§ Frequency of SSRs: per motif and number of repeats §§§# print OUT "Frequency of identified SSR motifs\n---------------------------- ------\n\nRepeats"; for ($i = $min_repeats;$i <= $max_repeats; $i++) {print OUT "\t$i"}; print OUT "\ttotal\n"; for ($i = 0; $i < scalar (@count_motifs); $i++) {my $typ = length ($count_motifs[$i]); print OUT $count_motifs[$i]; for ($j = $min_repeats; $j <= $max_repeats; $j++) {if ($j < $typrep{$typ}) {print OUT "\t-";next}; if ($count_motifs[$i]->{$j}) {print OUT "\t$count_motifs[$i]->{$j}"} else {print OUT "\t"}; }; print OUT "\t$count_motifs{$count_motifs[$i]}\n";}; print OUT "\n"; #§§§ Frequency of SSRs: summarizing redundant and reverse motifs §§§# # Eliminates %count_motifs ! print OUT "Frequency of classified repeat types (considering sequence complementary)\n----------------------------------------------------------- --------------\n\nRepeats"; my (%red_rev,@red_rev); # groups for ($i = 0; $i < scalar (@count_motifs); $i++) { next if ($count_motifs{$count_motifs[$i]} eq 'X'); my (%group,@group,$red_rev); # store redundant/reverse motifs 74 my $reverse_motif = $actual_motif = $actual_motif_a = $count_motifs[$i]; $reverse_motif =~ tr/ACGT/TGCA/; my $reverse_motif_a = $reverse_motif; for ($j = 0; $j < length ($count_motifs[$i]); $j++) {if ($count_motifs{$actual_motif}) {$group{$actual_motif} = "1"; $count_motifs{$actual_motif}='X'}; if ($count_motifs{$reverse_motif}) {$group{$reverse_motif} = "1"; $count_motifs{$reverse_motif}='X'}; $actual_motif =~ s/(.)(.*)/$2$1/; $reverse_motif =~ s/(.)(.*)/$2$1/; $actual_motif_a = $actual_motif if ($actual_motif lt $actual_motif_a); $reverse_motif_a = $reverse_motif if ($reverse_motif lt $reverse_motif_a) }; if ($actual_motif_a lt $reverse_motif_a) {$red_rev = "$actual_motif_a/$reverse_motif_a"} else {$red_rev = "$reverse_motif_a/$actual_motif_a"}; # group name $red_rev{$red_rev}++; @group = keys %group; for ($j = 0; $j < scalar (@group); $j++) {for ($k = $min_repeats; $k <= $max_repeats; $k++) { if ($group[$j]->{$k}) {$red_rev->{"total"} += $group[$j]- >{$k};$red_rev->{$k} += $group[$j]->{$k}} } }; for ($i = $min_repeats; $i <= $max_repeats; $i++) {print OUT "\t$i"}; print OUT "\ttotal\n"; @red_rev = sort {length ($a) length ($b) || $a cmp $b } keys %red_rev; for ($i = 0; $i < scalar (@red_rev); $i++) {my $typ = (length ($red_rev[$i])-1)/2; print OUT $red_rev[$i]; for ($j = $min_repeats; $j <= $max_repeats; $j++) {if ($j < $typrep{$typ}) {print OUT "\t-";next}; if ($red_rev[$i]->{$j}) {print OUT "\t",$red_rev[$i]->{$j}} else {print OUT "\t"}}; print OUT "\t",$red_rev[$i]->{"total"},"\n"; }; 3. Mã chƣơng trình ssrfinder_1_0  1_ssr_repeat_finder.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info $runtype = 1; # 1 = genbank fasta, 0 = local (fasta header differences) 75 # no need to change anything below this point $seqcount = 0; $ssr_count = 0; @ssr_label = qw(a b c d e f g h i j k l m n o p q r s t u v w x y z aa ab ac ad ae af ag ah ai aj ak al am an ao ap aq ar as at au av aw ax ay az); # set the flanking length here ( in bases ) $flank_length = 150; # i.e. 150 bp # set the minimum length of the repeat $min_pattern_length = 12; # i.e. 12 bp # open the input sequence file - fasta format 3rd field of header is the accession number open (SEQFILE, "../$datename/sequence$datename.txt") || die "file not found: $!"; # open the output file for the sequence ids open (IDFILE, ">>../$datename/new_ids$datename.txt") || die "couldn't create file"; # open the ouput file for the ssr results open (SSROUTFILE, ">>../$datename/ssrout$datename.txt") || die "couldn't create file"; #open the output file for the overall output for loading into labdb open (LABDBTXT, ">>../$datename/labdbout$datename.txt") || die "couldn't create file"; # read in file of previously checked IDs $CheckedIDs = `cat ../$datename/CheckedIDs.txt`; # parse the sequence file and process if ($runtype == 1) { while (defined ($line = )) { chomp $line; if ($line =~ /^>/) { if (defined ($Seq)) { # check if genbank id has been done before if ($CheckedIDs !~ /$SeqHead[3]/) { &SSRSearch; print "new seq\n"; } else { print "old seq\n";} } $seqcount++; print $seqcount; $HeadLine = $line; undef $Seq; @SeqHead = split(/\|/,$HeadLine); print IDFILE "$SeqHead[3]\t$SeqHead[4]\n"; 76 } else { $Seq = "$Seq" . "$line"; } if (defined ($Seq)) { # check if genbank id has been done before if ($CheckedIDs !~ /$SeqHead[3]/) { &SSRSearch; print "new seq\n"; } else { print "old seq\n"; } } } elsif ($runtype == 2) { while (defined ($line = )) { chomp ($line); @LineIn = split(/\t/, $line); $SeqHead[3] = $LineIn[0]; print $SeqHead[3]; $Seq = $LineIn[1]; if (defined ($Seq)) { &SSRSearch;} $seqcount++; print $seqcount; if ( ($seqcount%5) == 0 ) { print "\n"; } else { print "\t"; } undef $Seq; @SeqHead = split(/\|/,$HeadLine) } } elsif ($runtype == 3) { while (defined ($line = )) { chomp ($line); if ($line =~ /^>/) { if (defined ($Seq)) { &SSRSearch; } $seqcount++; print $seqcount; if ( ($seqcount%5) == 0 ) { print "\n"; } else { print "\t";} undef $Seq; $HeadLine = $line; @templine1 = split(/>/,$HeadLine); @templine2 = split(/\./, $templine1[1]); 77 $SeqHead[3] = $templine2[0]."_".$templine2[3]; print IDFILE "$SeqHead[3]\n"; } else { $Seq = "$Seq" . "$line"; } if (defined ($Seq)) { &SSRSearch; } close (SEQFILE); close (IDFILE); close (SSROUTFILE); close (LABDBTXT); print "Number of sequences in input file = $seqcount \n"; print "Number of Repeats found = $ssr_count \n"; exit 0; # subroutines sub SSRSearch() { print "*"; $suffix = -1; while ( $Seq =~ /(([ATGC]{2,})\2{3,})/gi ) { print "+"; $fullmatch = $1; $minmatch = $2; # minimum matches: di-nt 6 repeats, tri-nt 4 repeats, tetra-nt 4 repeats # check that repeat is greater than minimum total length of match # and that it is not a single nt repeat - AAAAAAA, TTTTTTTT, etc. if ( (( $length_sub = length($fullmatch)) >= $min_pattern_length ) && !( $fullmatch =~ /([ATGC])\1{8,}/ ) ) { print "-"; $ssr_count++; $suffix++; $Accession = "$SeqHead[3]" . "$ssr_label[$suffix]"; print SSROUTFILE "$Accession\t$fullmatch\t$minmatch\t"; print LABDBTXT "$SeqHead[3]\t$Accession\t$fullmatch\t$minmatch\t"; $pos2 = index($Seq,$fullmatch); $pos3 = $pos2 + $length_sub; if ($pos2 < $flank_length) { $pos1 = 0; } else { $pos1 = $pos2 - $flank_length; } if ( ( ( $max = length($Seq)) - $pos3 ) < $flank_length ) { $pos4 = $max; } else { $pos4 = $pos3 + $flank_length; } 78 $seqleft = substr($Seq,$pos1,$pos2-$pos1); $seqcenter = substr($Seq,$pos2,$pos3-$pos2); $seqright = substr($Seq,$pos3,$pos4-$pos3); print SSROUTFILE $seqleft, "[", $seqcenter,"]",$seqright,"\n"; print LABDBTXT $pos2,",",$length_sub,"\t",$seqleft, "[", $seqcenter,"]",$seqright,"\n"; } print "\n";}  2_ssr_primer_designer.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info $primer3app = 'd:\\detai\\ssrfinder_1_0\\primer3'; # full path to the primer3 command # no need to change anything below this point # open the ouput file with the ssr results open (SSRINFILE, "../$datename/ssrout$datename.txt") || die "couldn't open file"; open (RAWPRIMER, ">../$datename/raw_primer3$datename.out") || die "couldn't open file: $!"; open (PRIMEROUT, ">../$datename/primer_results$datename.txt") || die "couldn't open file: $!"; $counter = 0; while (defined($line = )) { chomp ($line); @columns = split(/\t/,$line); @seq = split(/\[|\]/,$columns[3]); $length_left = length($seq[0]); $length_mask = length($seq[1]); $length_right = length($seq[2]); $t1 = $length_left+1; $t2 = $length_mask; open (PRIMERIN, ">../$datename/primerin.txt") || die "couldn't open file: $!"; print PRIMERIN "PRIMER_SEQUENCE_ID=",$columns[1],"\n"; print PRIMERIN "SEQUENCE=",$seq[0],$seq[1],$seq[2],"\n"; print PRIMERIN "TARGET=", $t1, ",", $t2, "\n"; print PRIMERIN "PRIMER_PRODUCT_SIZE_RANGE=80-160 80-240 80-300\n"; print PRIMERIN "PRIMER_OPT_SIZE=24\n"; print PRIMERIN "PRIMER_MIN_SIZE=20\n"; print PRIMERIN "PRIMER_MAX_SIZE=28\n"; print PRIMERIN "PRIMER_OPT_TM=63\n"; 79 print PRIMERIN "PRIMER_MIN_TM=60\n"; print PRIMERIN "PRIMER_MAX_TM=65\n"; print PRIMERIN "PRIMER_MAX_DIFF_TM=1\n"; print PRIMERIN "=\n"; close (PRIMERIN); # $primer3 = `../devel/primer/primer3_0_9_test/src/primer3_core < ../$datename/primerin.txt`; $primer3 = `$primer3app < ../$datename/primerin.txt`; print RAWPRIMER "######## $columns[1] #########\n"; print RAWPRIMER $primer3, "\n"; @prime_out = split(/\n/, $primer3); foreach $i (0..$#prime_out) { ($varname,$varvalue) = split(/=/, $prime_out[$i]); $primehash{$varname} = $varvalue;} if ($primehash{'PRIMER_LEFT_SEQUENCE'}) { $counter++; print PRIMEROUT "$columns[0]\t$columns[1]\t$columns[2]"; print PRIMEROUT "\t$seq[0]$seq[1]$seq[2]"; print PRIMEROUT "\t", $primehash{'PRIMER_LEFT_SEQUENCE'}; print PRIMEROUT "\t", $primehash{'PRIMER_LEFT_TM'}; print PRIMEROUT "\t", $primehash{'PRIMER_RIGHT_SEQUENCE'}; print PRIMEROUT "\t", $primehash{'PRIMER_RIGHT_TM'}; print PRIMEROUT "\t", $primehash{'PRIMER_PRODUCT_SIZE'}; print PRIMEROUT "\n"; } undef %primehash;} close (PRIMEROUT); close (SSRINFILE); close (RAWPRIMER); print "repeats primed: ", $counter, "\n"; exit 0;  3_ssr_primer_rep_check.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info # no need to change anything below this point open (OUTFILE, ">rescreened$datename.txt"); $good = 0; $bad = 0; $temp = `type primer_results$datename.txt`; @temp1 = split(/\n/, $temp); foreach $i (0..$#temp1) { 80 @temp2 = split(/\t/, $temp1[$i]); print $temp2[0],"\t",$temp2[1],"\t",$temp2[2],"\t",$temp2[3]; print "\t",$temp2[4],"\t",$temp2[5],"\t",$temp2[6],"\t",$temp2[7],"\t",$temp2[8], "\t"; if (( $temp2[4] =~ /(([ATGC]{2,3})\2{3,})/gi ) || ( $temp2[6] =~ /(([ATGC]{2,3})\2{3,})/gi )) { print "bad\n"; $bad++; } else { print "good\n"; $good++; print OUTFILE $temp2[0],"\t",$temp2[1],"\t",$temp2[2],"\t",$temp2[3],"\t"; print OUTFILE $temp2[4],"\t",$temp2[5],"\t",$temp2[6],"\t",$temp2[7],"\t",$temp2[8],"\n" } print "good: $good\nbad: $bad\n"; close (OUTFILE); exit(0);  4_ssr_primer_blast.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info $blastapp = 'd:\\detai\\ssrfinder_1_0\\blastall'; # full path to the blastall command $blastdbdir = 'd:\\detai\\ssrfinder_1_0\\db'; #full path to the blast database directory $blastdbname = 'AllPrimers.nt'; # name of the blast database to use $formatdbapp = 'd:\\detai\\ssrfinder_1_0\\formatdb'; # full path to the formatdb command # no need to change anything below this point open (PRIMERSIN, "../$datename/rescreened$datename.txt") || die; open (BLASTOUTFILE, ">../$datename/blastout$datename.txt") || die; open (FULLBLASTOUT, ">../$datename/fullblastoutput$datename.txt") || die; $counter = 0; $b0hit = 0; while (defined($line=)) { chomp ($line); @columns = split(/\t/,$line); $Accession = $columns[0]; $sequence = $columns[3]; 81 $forward = $columns[4]; $reverse = $columns[6]; $counter++; print $counter, "\n"; print BLASTOUTFILE $Accession, "\t", $columns[1], "\t", $columns[2]; print BLASTOUTFILE "\t", $sequence, "\t", $forward, "\t", $columns[5], "\t"; print BLASTOUTFILE $reverse, "\t", $columns[7], "\t", $columns[8], "\t"; &BlastIt; print "Blasted Sequences= ",$counter,"\n"; print "Blast NON-Hits= ",$b0hit,"\n"; close (PRIMERSIN); close (BLASTOUTFILE); close (FULLBLASTOUT); exit 0; sub BlastIt() { open (TMPSEQFILE, ">../$datename/repeats/$Accession.fasta"); print TMPSEQFILE "> ",$Accession,"\n",$sequence,"\n"; close (TMPSEQFILE); # $blastout = `blastall -p blastn -d db/AllPrimers.nt -e 0.01 - i ../$datename/repeats/$Accession.fasta`; $blastout = `$blastapp -p blastn -d $blastdbdir/$blastdbname -e 0.01 -i ../$datename/repeats/$Accession.fasta`; print FULLBLASTOUT "XXXXX\t",$Accession,"\tXXXXXXXXXXXXXXXXXXXXXXXXX\n\n"; print FULLBLASTOUT $blastout; &ParseBlast; sub DBup() { # open (BLASTDBFASTA, ">>../blast/db/AllPrimers.nt"); open (BLASTDBFASTA, ">>$blastdbdir/$blastdbname"); print BLASTDBFASTA "> ",$Accession,"_f\n",$forward,"\n"; print BLASTDBFASTA "> ",$Accession,"_r\n",$reverse,"\n"; close (BLASTDBFASTA); # $formatdbstatus = system("formatdb -i db/AllPrimers.nt -p F -o T"); $formatdbstatus = system("$formatdbapp -i $blastdbdir/$blastdbname -p F -o T");} sub ParseBlast() { # split the blast output into records (splits at blank lines) @blastrecs = split(/\n\n/, $blastout); # the record that contains the hit results SHOULD be in record 6 - WATCH OUT if ($blastrecs[5] =~ /Sequence/) { 82 @blastmatch = split(/\n/, $blastrecs[6]); print BLASTOUTFILE $#blastmatch-1, "\t"; for $entry (0 .. $#blastmatch) { ($name,$score, $E) = split(/[ ]+ /, $blastmatch[$entry]); print $name,"\t",$score,"\t",$E,"\n"; print BLASTOUTFILE "\t",$name,"\t",$score,"\t",$E;} } else { print BLASTOUTFILE "0\t"; $b0hit++; print "No Hits found\n"; print BLASTOUTFILE "\tnone\tnull\tnull"; &DBup; } print BLASTOUTFILE "\n";}  5_ssr_order_filter.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info # no need to change anything below this point open (OUTFILE, ">filter$datename.txt"); $infile = `cat blastout$datename.txt`; $count = 0; @temp1 = split(/\n/, $infile); foreach $f (0..$#temp1) { @temp2 = split(/\t/, $temp1[$f]); if ($temp2[11] eq 'none') { print OUTFILE $temp2[0],"\t",$temp2[1],"\t"; print OUTFILE $temp2[2],"\t",$temp2[3],"\t"; print OUTFILE $temp2[4],"\t",$temp2[5],"\t"; print OUTFILE $temp2[6],"\t",$temp2[7],"\t"; print OUTFILE $temp2[8],"\n"; print $temp2[0],"\t",$temp2[10],"\t",$temp2[11],"\n"; $count++;} print "non-hits: ",$count,"\n"; close (OUTFILE); exit (0);  6_ssr_order_formatter.pl #!/usr/bin/perl –w # change these parameters for each run!!!!!!!! $datename = '20060715'; # date for directory name and datafile info 83 # no need to change anything below this point open (OUTFILE, ">order$datename.txt"); $infile = `cat blastout$datename.txt`; $count = 0; @temp1 = split(/\n/, $infile); foreach $f (0..$#temp1) { @temp2 = split(/\t/, $temp1[$f]); if ($temp2[11] eq 'none') { print OUTFILE $temp2[0],"\t"; print OUTFILE $temp2[4],"\t",$temp2[5],"\t"; print OUTFILE $temp2[6],"\t",$temp2[7],"\t"; print OUTFILE $temp2[8],"\n"; print $temp2[0],"\t",$temp2[10],"\t",$temp2[11],"\n"; $count++; }} print "non-hits: ",$count,"\n"; close (OUTFILE); exit (0);