-
Notifications
You must be signed in to change notification settings - Fork 0
/
FragmentBuild.m
50 lines (38 loc) · 1.23 KB
/
FragmentBuild.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
clear all
close all
FragmentLength = 1000;
UnitRegionSize = 100000;
[Head Seq] = fastaread('chromFA/chr22.fa');
fid = fopen(['DomainStartsAndEnds.txt'], 'r');
tline = fgetl(fid);
i = 0;
while ischar(tline)
i = i+1;
C = strsplit(tline);
DomainStart(i) = str2num(C{2});
DomainEnd(i) = str2num(C{3});
tline = fgetl(fid);
end
fclose(fid);
for i = 1:length(DomainStart)
DomainCenter(i) = (DomainStart(i)+DomainEnd(i))/2;
if DomainCenter(i)-UnitRegionSize/2 >= DomainStart(i)
RegionStart(i) = DomainCenter(i)-UnitRegionSize/2;
else
RegionStart(i) = DomainStart(i);
end
if DomainCenter(i)+UnitRegionSize/2 <= DomainEnd(i)
RegionEnd(i) = DomainCenter(i)+UnitRegionSize/2;
else
RegionEnd(i) = DomainEnd(i);
end
Domain{i} = Seq(round(RegionStart(i)):round(RegionEnd(i)));
end
for i = 1:length(DomainStart)
for j = 1:floor(length(Domain{i})/FragmentLength)
data(j).Sequence = Domain{i}((j-1)*FragmentLength+1:j*FragmentLength);
data(j).Header = ['Domain' num2str(i) 'Frag' num2str(j)];
end
fastawrite(['1kb_fragments_sequence' num2str(i) '.fasta'], data);
clear data
end