Final Round. Problem 1 Gene Expression
https://stepik.org/lesson/32224/step/2?unit=19659
Problem 1
Problem 1
stepik.org
하나의 genome에서 유래한 expression level을 구하는 문제이다.
N개의 genome과 그것들로부터 유래된 rna 조각들을 M개를 주었을때 하나의 genome으로 부터 유래한 rna의 갯수를 세어야 한다. 즉 2개의 genome과 겹치는 rna들을 제외해야 한다.
간단하게 생각해보면 다음과 같은 코드로 가능하다.
import sys
from functools import reduce
import bisect
input = []
# with open('F1') as f:
# for line in f.readlines():
# input.append(line.strip())
for line in sys.stdin:
input.append(line.strip())
def read_to_intervals(genome, idx):
intervals = []
for i in range(int(len(genome)/2)):
intervals.append((genome[i*2], genome[i*2+1], idx))
return intervals
data = list(map(lambda x: list(map(lambda y: int(y), x.split(' '))), input))
n, m = data[0]
genomes = data[1: 1+n]
reads = data[-m:]
genomes = reduce(lambda a, b: a + b, (map(lambda genome: read_to_intervals(genome[1], genome[0]), enumerate(genomes))))
reads = list(map(lambda genome: read_to_intervals(genome[1], genome[0]), enumerate(reads)))
start_sorted = sorted(genomes, key=lambda x: x[0])
sorted_starts = list(map(lambda x: x[0], start_sorted))
result = [0] * n
for i, read in enumerate(reads):
intersect_genes = set()
for interval in read:
read_start, read_end, idx = interval
start_idx = bisect.bisect(sorted_starts, read_end)
candidates = start_sorted[0:start_idx]
end_sorted = sorted(candidates, key=lambda x: x[1])
sorted_ends = list(map(lambda x: x[1], end_sorted))
end_idx = bisect.bisect_left(sorted_ends, read_start)
intersects = end_sorted[end_idx:]
intersect_genes = intersect_genes.union(set(map(lambda intersect: intersect[2], intersects)))
if len(intersect_genes) > 1:
break
if len(intersect_genes) == 1:
result[intersect_genes.pop()] += 1
print('\n'.join(map(lambda x: str(x), result)))genome들의 read와 rna들의 read를 각각 정렬한후 rna의 raad들을 순회하며 해당 read와 겹치는 genome들의 index를 set에 저장하는 것이다. 해당 rna와 겹치는 dna 탐색을 조금 더 빠르게 하기위해 bisect, 2진탐색을 이용하였다.
rna의 read와 겹치는 판별은 다음과 같이 한다.
1. genome의 read들을 start index를 기준으로 정렬한다.
2. rna read의 end index보다 작은 start index 값을 가지는 genome read들을 선별해낸다.
3. 선별해낸 genome read들을 end index 기준으로 정렬한다.
4. rna read의 start index보다 큰 end index를 가지는 genome read들을 선별해낸다.
결과는 다음과 같이 산출한다.
1. 위 과정에서 선별해는 genome read들의 genome idx값을(0~N-1) 추출해낸다.
2. 해당 genome idx의 값들의 set을 rna idx를 key로하는 map에 저장한다.
3. rna idx를 key로하는 map을 순회하며 value에 있는 genome의 expression level에 1을 더한다. 만약 value에 있는 set의 길이가 1보다 크다면 건너뛴다.
위 과정의 시간복잡도를 생각해보자. 우선 genome read의 길이가 i이면 genome은 n개가 있으니 총 read는 in이고 정렬하는데 inlog(in)이 필요하다. 그리고 길이가 k인 rna read들을 순회하는데는 km의 필요하다. rna read의 end index보다 작은 start index 값을 가지는 genome read들을 선별해내는데 log(in)이 필요하며 이를 정렬하는데 또 inlog(in)이 필요하다. 이를 다 합해보면 O(in(log(in) + kmlog(km)*(log(in) + inlog(in)))이라고 볼 수 있다. 이는 당연히 easy에서도 타임아웃이 난다.
우리는 좀더 시간복잡도를 줄여야 한다. 그리고 며칠동안 고민한 끝에 모든 read들을 genome과 rna 구분없이 모아서 정렬 한 후 앞에서부터 하나씩 순회하면서 현재 열려있는 rna에 현재 열려있는 genome idx들을 저장하면 된다는 생각을 하게 되었다.
전체코드는 다음과 같다.
import sys
input = []
# with open('F1') as f:
# for line in f.readlines():
# input.append(line.strip())
for line in sys.stdin:
input.append(line.strip())
def read_to_intervals(genome, idx):
intervals = []
for i in range(int(len(genome)/2)):
intervals.append((genome[i*2], genome[i*2+1], idx))
return intervals
data = list(map(lambda x: list(map(lambda y: int(y), x.split(' '))), input))
n, m = data[0]
genomes = data[1: 1+n]
reads = data[-m:]
flatten_genomes = []
for idx, genome in enumerate(genomes):
for s_e, pos in enumerate(genome):
flatten_genomes.append([pos, idx, 'S' if s_e % 2 == 0 else 'E', 'G'])
flatten_reads = []
for idx, read in enumerate(reads):
for s_e, pos in enumerate(read):
flatten_reads.append([pos, idx, 'S' if s_e % 2 == 0 else 'E', 'R'])
gen_reads = sorted(flatten_genomes + flatten_reads, key= lambda x: x[2], reverse=True)
gen_reads = sorted(gen_reads, key= lambda x: x[0])
prev_genes_len = 0
cur_gens = set()
cur_reads = set()
read_gen_map = [None for i in range(m)]
expressions = [0 for i in range(n)]
for idx, read in enumerate(gen_reads):
prev_genes_len = len(cur_gens)
if read[3] == 'G':
if read[2] == 'S':
cur_gens.add(read[1])
else:
cur_gens.remove(read[1])
else:
if read[2] == 'S' and read_gen_map[read[1]] != -1:
cur_reads.add(read[1])
elif read[2] == 'E':
if read[1] in cur_reads:
cur_reads.remove(read[1])
if len(cur_gens) > 1:
for cur_read in cur_reads:
read_gen_map[cur_read] = -1
cur_reads.clear()
elif len(cur_gens) == 1:
new_cur_reads = set()
if prev_genes_len == 1:
cur_read = read[1]
if read_gen_map[cur_read] is None:
read_gen_map[cur_read] = list(cur_gens)[0]
new_cur_reads.add(cur_read)
elif read_gen_map[cur_read] != list(cur_gens)[0]:
read_gen_map[cur_read] = -1
else:
new_cur_reads.add(cur_read)
else:
for cur_read in cur_reads:
if read_gen_map[cur_read] is None:
read_gen_map[cur_read] = list(cur_gens)[0]
new_cur_reads.add(cur_read)
elif read_gen_map[cur_read] != list(cur_gens)[0]:
read_gen_map[cur_read] = -1
else:
new_cur_reads.add(cur_read)
cur_reads = new_cur_reads
for gen_set in read_gen_map:
if gen_set is not None and gen_set != -1:
expressions[gen_set] += 1
print('\n'.join(map(lambda x: str(x), expressions)))차례차례 보도록 하자.
import sys
input = []
# with open('F1') as f:
# for line in f.readlines():
# input.append(line.strip())
for line in sys.stdin:
input.append(line.strip())
def read_to_intervals(genome, idx):
intervals = []
for i in range(int(len(genome)/2)):
intervals.append((genome[i*2], genome[i*2+1], idx))
return intervals
data = list(map(lambda x: list(map(lambda y: int(y), x.split(' '))), input))
n, m = data[0]
genomes = data[1: 1+n]
reads = data[-m:]
flatten_genomes = []
for idx, genome in enumerate(genomes):
for s_e, pos in enumerate(genome):
flatten_genomes.append([pos, idx, 'S' if s_e % 2 == 0 else 'E', 'G'])
flatten_reads = []
for idx, read in enumerate(reads):
for s_e, pos in enumerate(read):
flatten_reads.append([pos, idx, 'S' if s_e % 2 == 0 else 'E', 'R'])
# flatten_genomes = reduce(lambda a,b: a+b, list(map(lambda idx_gene: list(map(lambda idx_read: [idx_read[1], idx_gene[0], 'S' if idx_read[0] % 2 == 0 else 'E', 'G'], enumerate(idx_gene[1]))), enumerate(genomes))))
# flatten_reads = reduce(lambda a,b: a+b, list(map(lambda idx_gene: list(map(lambda idx_read: [idx_read[1], idx_gene[0], 'S' if idx_read[0] % 2 == 0 else 'E', 'R'], enumerate(idx_gene[1]))), enumerate(reads))))
gen_reads = sorted(flatten_genomes + flatten_reads, key= lambda x: x[2], reverse=True)
gen_reads = sorted(gen_reads, key= lambda x: x[0])
위 코드는 단순하게 genome와 rna의 read들을 모두 하나의 list로 만들어 정렬하는 것이다. 다만 genome은 'G', rna는 'R', 시작은 'S', 끝은 'E'로 표시하여 [pos, genome idx or rna idx, 'S' or 'E', 'G' or 'R']로 만든다. 이때 for문으로 순회할때 E가 먼저 나오게 되면 pos가 start와 end가 겹치는 경우가 발생하는 경우가 빠지게 되므로 S가 E보다 먼저 나와야 read가 닫히기 전에 다른 read가 시작되어 해당 경우가 정확하게 count 될수 있다. 그래서 마지막에 sort를 한번은 'E', 'S'로 정렬해 준 후 pos로 정렬하는 것이다. python의 sort함수는 안정정렬이므로 위와같이 하면 된다.
prev_genes_len = 0
cur_gens = set()
cur_reads = set()
read_gen_map = [None for i in range(m)]
expressions = [0 for i in range(n)]prev_genes_len이 있는 이유는 뒤에 살펴보도록 한다.
cur_gens와 cur_reads는 현재 열려있는 genome과 rna의 set이다. 'S'면 set에 추가하고 'E'면 set에서 제거한다.
read_gen_map은 해당 rna의 read가 어느 genome에서 왔는지 genome의 idx를 저장하는 배열이다. i의 rna가 j의 dna로 부터 유래했다면 read_gen_map[i] = j가 된다. 만약 두개 이상의 genome에서 유래했다면 -1을 저장한다.
expressions는 i번째 genome의 expression level을 저장하는 배열이다.
if read[3] == 'G':
if read[2] == 'S':
cur_gens.add(read[1])
else:
cur_gens.remove(read[1])
else:
if read[2] == 'S' and read_gen_map[read[1]] != -1:
cur_reads.add(read[1])
elif read[2] == 'E':
if read[1] in cur_reads:
cur_reads.remove(read[1])'G'이고 'S'면 cur_gens에 추가. 'E'면 cur_gens에서 제거
'R'이고 'S'먄 cur_reads에 추가. 'E'먄 cur_gens에서 제거. 다만 이미 여러개의 genome과 겹친다고 판명난 rna는 괜히 또 추가하지 않는다.
if len(cur_gens) > 1:
for cur_read in cur_reads:
read_gen_map[cur_read] = -1
cur_reads.clear()현재 열려있는 genome이 1개보다 많으면 현재 열려있는 rna의 genome 매핑에 -1을 넣는다. 그리고 현재 reads들은 전부 다 제거한다. 괜히 해당 rna의 read들이 닫힐때까지 저장하지 않아도 되기 때문이다.
elif len(cur_gens) == 1:
new_cur_reads = set()
if prev_genes_len == 1:
cur_read = read[1]
if read_gen_map[cur_read] is None:
read_gen_map[cur_read] = list(cur_gens)[0]
new_cur_reads.add(cur_read)
elif read_gen_map[cur_read] != list(cur_gens)[0]:
read_gen_map[cur_read] = -1
else:
new_cur_reads.add(cur_read)
else:
for cur_read in cur_reads:
if read_gen_map[cur_read] is None:
read_gen_map[cur_read] = list(cur_gens)[0]
new_cur_reads.add(cur_read)
elif read_gen_map[cur_read] != list(cur_gens)[0]:
read_gen_map[cur_read] = -1
else:
new_cur_reads.add(cur_read)
cur_reads = new_cur_reads
만일 현재 열려있는 genome이 1개라면 여러가지 경우가 있다.
1. 현재 열려있는 rna가 아직 아무런 genome에 매핑되지 않은 경우 -> read_gen_map의 현재 rna에 현재 genome 매핑 추가
2. 현재 열려있는 rna에 이미 다른 genome 매핑이 돼있는 경우 -> read_gen_map의 현재 rna에 -1 매핑
3. 현재 열려있는 rna에 같은 genome이 매핑이 돼있는 경우 -> pass
new_cur_reads는 -1이 된 rna들을 제거하기 위한 것이다.
위 설명은 prev_genes_len의 역할이 있지 않다. 실제로 위 로직으로 hard에 제출한경우 10번째 케이스에서 타임아웃이 나게된다. 이유는 현재 열린 genome이 1개인 상태로(cur_gens의 길이가 1인 상태로) 엄청나게 많은 rna read들이 동시에 열리게 되면서 cur_reads에 수많은 rna idx가 들어있게되어 for문을 계속해서 순회하기 때문이다. 만일 read를 처리하기 전과 후과 똑같이 cur_gens의 길이가 1이라면 현재의 read는 rna의 read라는걸 생각할 수 있다. genome의 read라면 길이가 0이 되거나 2가 돼야하기 때문이다. 그러므로 현재 read의 rna만 처리하면 되므로 for문을 돌면서 열려있는 모든 rna idx를 순회할 필요가 없어진다.