comparison rDiff/src/tools/read_utils/splicingsequence.m @ 0:0f80a5141704

version 0.3 uploaded
author vipints
date Thu, 14 Feb 2013 23:38:36 -0500
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:0f80a5141704
1 function [SPLICINGEVENTS, SEQUENCE, EXONSEQUENCE, IDENTIFICATIONLENGTH] = splicingsequence(GENE)
2 % This function generates all sequence of all splicesites
3 % SPLICINGEVENTS, a SEQUENCE marking which elements in
4 % SPLICINGEVENTS are introns and exons and a sequence for each
5 % transcript which indicates which Exons are in cluded in this transcript
6 % SEQUENCE contains indices of SPLICINGEVENTS, containing 1, if the position right to the idx in the transcript
7 % are exonic, 0 otherwise
8
9 EXONS = GENE.exons;
10 START = GENE.start;
11 STOP = GENE.stop;
12
13 NB_OF_TRANSCR = size(EXONS,2);
14
15 SPLICINGEVENTS = [];
16 for i = 1:NB_OF_TRANSCR
17 SPLICINGEVENTS = [SPLICINGEVENTS, EXONS{i}(:,1)', EXONS{i}(:,2)' ];
18 end
19 SPLICINGEVENTS = SPLICINGEVENTS - START + 1;
20 SPLICINGEVENTS = unique(SPLICINGEVENTS);
21
22 EXONSEQUENCE = zeros(NB_OF_TRANSCR, length(SPLICINGEVENTS) - 1);
23
24 for i = 1:NB_OF_TRANSCR
25 %%% for every exon in transcript i
26 for j = 1:size(EXONS{i}, 1)
27 POS_START = find(SPLICINGEVENTS == EXONS{i}(j,1) - START + 1, 1, 'first');
28 POS_STOP = find(SPLICINGEVENTS < EXONS{i}(j,2) - START + 1, 1, 'last');
29 %%% set splicing events active in transcript i to 1
30 EXONSEQUENCE(i, POS_START:POS_STOP) = 1;
31 end
32 end
33
34 %%% merge active splicing events of all transcripts
35 SEQUENCE = max(EXONSEQUENCE, [], 1);
36 IDENTIFICATIONLENGTH = [];
37