#python #bioinformatics
#python #биоинформатика
Вопрос:
У меня есть несколько файлов образцов (> 20), которые выглядят как:
chr startpos endpos
1 14930 818094
1 818161 31595422
2 35593931 35865807
2 35868158 104785784
И я хотел бы вывести области, которые являются общими для образцов. Например. если образец 1 имеет:
1 14900 818000
пример 2:
1 15000 605000
пример 3:
1 25000 705000
Я хотел бы вывести:
1 25000 605000
Я также хотел бы включить правило большинства, такое, что, например, если 10 из 20 образцов имеют минимальный регион -> вывести регион. Т.е. я хотел бы гибко указать, сколько образцов должны иметь регион, чтобы он был напечатан на выходе.
У кого-нибудь есть решение на python для этого?
Комментарии:
1. На самом деле этот вопрос касается не Unix / Linux, а программирования (кодирования, алгоритмов), поэтому он больше подходит для Stack Overflow, а не для этого сайта.
2. Также обратите внимание, что сотрудники Stack Exchange обычно не хотят выполнять вашу работу / домашнее задание за вас. Здесь есть добровольцы, которые рады помочь, но вам нужно показать, что вы тоже прилагаете усилия. Поэтому попробуйте решить это самостоятельно и, когда вы зайдете в тупик, задайте конкретный вопрос о том, что происходит неожиданно. У вас больше шансов получить полезные ответы (и научиться!) таким образом.
3. Да, извините, моя ошибка
4. поиск в Google по «bedtools» — bedtools intersect может быть хорошим местом для начала
Ответ №1:
Не уверен, относится ли этот вопрос к Unix amp; Linux stackexchange. Это больше похоже на общий вопрос программирования.
Тем не менее, я бы посоветовал вам изучить возможность использования pandas .
Вы можете импортировать файл образца в виде фрейма данных, указав следующие обозначения табуляции:
import pandas as pd
df = pd.read_csv('/tmp/samplefile.csv',sep='t')
Если вы знаете, что начальные значения всегда будут меньше конечных значений, вы могли бы найти искомый результат, взяв максимум df['startpos'] и минимум df['endpos'] .
Комментарии:
1. Спасибо за начало, я попробую
Ответ №2:
Существуют инструменты, предназначенные специально для работы с данными формата BED, изучение которых вам было бы полезно. bedtools, вероятно, самый распространенный и простой в использовании, и у него есть оболочка python, доступная, если вам абсолютно необходимо использовать python.
Инструмент с несколькими разделениями — это то, что вам, вероятно, нужно, и он должен быть довольно простым в использовании:
Пример использования:
== Input files: ==
$ cat a.bed
chr1 6 12
chr1 10 20
chr1 22 27
chr1 24 30
$ cat b.bed
chr1 12 32
chr1 14 30
$ cat c.bed
chr1 8 15
chr1 10 14
chr1 32 34
$ cat sizes.txt
chr1 5000
== Multi-intersect the files: ==
$ multiIntersectBed -i a.bed b.bed c.bed
chr1 6 8 1 1 1 0 0
chr1 8 12 2 1,3 1 0 1
chr1 12 15 3 1,2,3 1 1 1
chr1 15 20 2 1,2 1 1 0
chr1 20 22 1 2 0 1 0
chr1 22 30 2 1,2 1 1 0
chr1 30 32 1 2 0 1 0
chr1 32 34 1 3 0 0 1
== Multi-intersect the files, with a header line (titles are the file names): ==
$ multiIntersectBed -header -i a.bed b.bed c.bed
chrom start end num list a.bed b.bed c.bed
chr1 6 8 1 1 1 0 0
chr1 8 12 2 1,3 1 0 1
chr1 12 15 3 1,2,3 1 1 1
chr1 15 20 2 1,2 1 1 0
chr1 20 22 1 2 0 1 0
chr1 22 30 2 1,2 1 1 0
chr1 30 32 1 2 0 1 0
chr1 32 34 1 3 0 0 1
== Multi-intersect the files, with a header line and custom names: ==
$ multiIntersectBed -header -i a.bed b.bed c.bed -names A B C
chrom start end num list A B C
chr1 6 8 1 A 1 0 0
chr1 8 12 2 A,C 1 0 1
chr1 12 15 3 A,B,C 1 1 1
chr1 15 20 2 A,B 1 1 0
chr1 20 22 1 B 0 1 0
chr1 22 30 2 A,B 1 1 0
chr1 30 32 1 B 0 1 0
chr1 32 34 1 C 0 0 1
== Multi-intersect the files, showing empty regions (note, requires -g): ==
$ multiIntersectBed -header -i a.bed b.bed c.bed -names A B C -empty -g sizes.txt
chrom start end num list A B C
chr1 0 6 0 none 0 0 0
chr1 6 8 1 A 1 0 0
chr1 8 12 2 A,C 1 0 1
chr1 12 15 3 A,B,C 1 1 1
chr1 15 20 2 A,B 1 1 0
chr1 20 22 1 B 0 1 0
chr1 22 30 2 A,B 1 1 0
chr1 30 32 1 B 0 1 0
chr1 32 34 1 C 0 0 1
chr1 34 5000 0 none 0 0 0
Затем вы можете выполнить фильтрацию по 4-му столбцу, чтобы выполнить указанное вами условие для областей.