2015深圳杯数学建模-DNA序列k-mer index问题
2015-05-05 16:19
281 查看
这是一个数学建模的问题,非常考察编程能力。
我们选择这个题目,就是因为i这个题目可以转化为模式串匹配问题。
由于这个题目中的模式串为k-mer,数量特别大,我们队用了多模式匹配的AC自动机算法找位置。
但是存储k-mer位置的索引表的数据量也很大,我们就用了哈希函数的直接定址法+链表的方式建立索引表。
虽然在海量数据下,建立索引的时间很长,但是一旦建立索引再去查找位置非常快。
题目如下:
这个问题来自 DNA序列的k-mer index问题。
给定一个DNA序列,这个系列只含有4个字母ATCG,如 S =“CTGTACTGTAT”。给定一个整数值k,从S的第一个位置开始,取一连续k个字母的短串,称之为k-mer(如k= 5,则此短串为CTGTA), 然后从S的第二个位置, 取另一k-mer(如k= 5,则此短串为TGTAC),这样直至S的末端,就得一个集合,包含全部k-mer 。 如对序列S来说,所有5-mer为
{CTGTA,TGTAC,GTACT,TACTG,ACTGT,TGTAT}
通常这些k-mer需一种数据索引方法,可被后面的操作快速访问。例如,对5-mer来说,当查询CTGTA,通过这种数据索引方法,可返回其在DNA序列S中的位置为{1,6}。
问题
现在以文件形式给定 100万个 DNA序列,序列编号为1-1000000,每个基因序列长度为100 。
(1)要求对给定k, 给出并实现一种数据索引方法,可返回任意一个k-mer所在的DNA序列编号和相应序列中出现的位置。每次建立索引,只需支持一个k值即可,不需要支持全部k值。
(2)要求索引一旦建立,查询速度尽量快,所用内存尽量小。
(3)给出建立索引所用的计算复杂度,和空间复杂度分析。
(4)给出使用索引查询的计算复杂度,和空间复杂度分析。
(5)假设内存限制为8G,分析所设计索引方法所能支持的最大k值和相应数据查询效率。
(6)按重要性由高到低排列,将依据以下几点,来评价索引方法性能
· 索引查询速度
· 索引内存使用
· 8G内存下,所能支持的k值范围
· 建立索引时
我们队高大上的代码:
这是一个数学建模的问题,非常考察编程能力。
我们选择这个题目,就是因为i这个题目可以转化为模式串匹配问题。
由于这个题目中的模式串为k-mer,数量特别大,我们队用了多模式匹配的AC自动机算法找位置。
但是存储k-mer位置的索引表的数据量也很大,我们就用了哈希函数的直接定址法+链表的方式建立索引表。
虽然在海量数据下,建立索引的时间很长,但是一旦建立索引再去查找位置非常快。
题目如下:
这个问题来自 DNA序列的k-mer index问题。
给定一个DNA序列,这个系列只含有4个字母ATCG,如 S =“CTGTACTGTAT”。给定一个整数值k,从S的第一个位置开始,取一连续k个字母的短串,称之为k-mer(如k= 5,则此短串为CTGTA), 然后从S的第二个位置, 取另一k-mer(如k= 5,则此短串为TGTAC),这样直至S的末端,就得一个集合,包含全部k-mer 。 如对序列S来说,所有5-mer为
{CTGTA,TGTAC,GTACT,TACTG,ACTGT,TGTAT}
通常这些k-mer需一种数据索引方法,可被后面的操作快速访问。例如,对5-mer来说,当查询CTGTA,通过这种数据索引方法,可返回其在DNA序列S中的位置为{1,6}。
问题
现在以文件形式给定 100万个 DNA序列,序列编号为1-1000000,每个基因序列长度为100 。
(1)要求对给定k, 给出并实现一种数据索引方法,可返回任意一个k-mer所在的DNA序列编号和相应序列中出现的位置。每次建立索引,只需支持一个k值即可,不需要支持全部k值。
(2)要求索引一旦建立,查询速度尽量快,所用内存尽量小。
(3)给出建立索引所用的计算复杂度,和空间复杂度分析。
(4)给出使用索引查询的计算复杂度,和空间复杂度分析。
(5)假设内存限制为8G,分析所设计索引方法所能支持的最大k值和相应数据查询效率。
(6)按重要性由高到低排列,将依据以下几点,来评价索引方法性能
· 索引查询速度
· 索引内存使用
· 8G内存下,所能支持的k值范围
· 建立索引时
我们队高大上的代码:
#include <algorithm> #include <iostream> #include <fstream> #include <cstring> #include <cstdlib> #include <cstdio> #include <cmath> #include <queue> #include <stack> #include <set> #include <map> using namespace std; const long long maxn = 1048580; const int INF = 0x3fffff; const double EPS = 1e-6; const double PI = 3.1415926; typedef long long ll; int dx[4] = {0,0,1,-1}; int dy[4] = {1,-1,0,0}; int k; typedef struct node //储值 { int add; struct node *next; } node,*Rlink; Rlink head[maxn]; void List_insert(int add,int data) //链表中插入地址 { if(head[add]==NULL) { head[add] = (node *)malloc(sizeof(node)); head[add]->add = data; head[add]->next = NULL; } else { node *p = head[add]; while(p->next) p=p->next; node *q = (node *)malloc(sizeof(node)); q->add = data; q->next = NULL; p->next = q; } } char str[1048580][101]; char lstr[110000000]; int hash(char *strr); struct Trie { int next[1048580][101],fail[1048580],end[1048580]; int root,L; int newnode() { for(int i = 0; i < 101; i++) next[L][i] = -1; end[L++] = -1; return L-1; } void init() { L = 0; root = newnode(); } void insert(char s[],int id) { int len = strlen(s); int now = root; for(int i = 0; i < len; i++) { if(next[now][s[i]] == -1) next[now][s[i]] = newnode(); now = next[now][s[i]]; } end[now] = id; } void build() { queue<int>Q; fail[root] = root; for(int i = 0; i < 101; i++) if(next[root][i] == -1) next[root][i] = root; else { fail[next[root][i]] = root; Q.push(next[root][i]); } while(!Q.empty()) { int now = Q.front(); Q.pop(); for(int i = 0; i < 101; i++) if(next[now][i] == -1) next[now][i]=next[fail[now]][i]; else { fail[next[now][i]]=next[fail[now]][i]; Q.push(next[now][i]); } } } int num[1048580];//记录个数 vector<int>post[1048580];//记录位置 void query(char buf[],int n) { for(int i = 0; i < n; i++) num[i] = 0,post[i].clear(); int len=strlen(buf); int now=root; for(int i=0; i<len; i++) { now=next[now][buf[i]]; int temp = now; while( temp != root ) { if(end[temp] != -1) { num[end[temp]]++; post[end[temp]].push_back(i); } temp = fail[temp]; } } for(int i = 0; i < n; i++) { if(num[i] > 0) { int add = hash(str[i]); //存储在数组中 for(int j=0; j<post[i].size(); j++) { List_insert(add,post[i][j]+1-k); // } } } } }; Trie ac; int hash(char *strr) { int len = strlen(strr); int add = 0; //地址 int base = 1; //进制 for(int i=len-1; i>=0; i--) { switch(strr[i]) { case 'A': add += 0*base; break; case 'T': add += 1*base; break; case 'G': add += 2*base; break; case 'C': add += 3*base; break; } base*=4; } return add; } void substr(char *st,char *ss,int i,int j) //i开始截取j个字符,扔进st中 { for(int k=0; k<j; k++) st[k] = ss[k+i]; st[j] = '\0'; } void Index(char *s) { int add = hash(s); node *p = head[add]; while(p) { printf("模式串所在DNA序列号为: %d,所在的下标为: %d\n",p->add/101+1,p->add%101+1); p=p->next; } } int main() { printf("请输入一个K值:\n"); scanf("%d",&k); map <string,int> map; //去重 string ss; char tmp[102]; ac.init(); int i=0; int cnt=0; //记录模式串个数 ifstream in("test.txt"); //文件读取 if(in) { while(getline(in,ss)) { ss.copy(tmp,100,0); tmp[100] = '*'; //将不同行的字符串分隔开,避免跨行匹配成功 tmp[101] = '\0'; strcat(lstr,tmp); for(int j=0; j<=100-k; j++) { string moshi = ss.substr(j,k); if(map[moshi]==1) continue; //已经插入过的模式串不会再插入 substr(str[cnt],tmp,j,k); ac.insert(str[cnt],cnt); cnt++; map[moshi] = 1; //插入标记 } i++; } } ac.build(); //ac自动机建立 ac.query(lstr,cnt); char l[101]; printf("请输入您要查询的k-mer:\n"); while(~scanf("%s",l)) Index(l); return 0; }
相关文章推荐
- 2015年“深圳杯”数学建模夏令营-B题:DNA序列的k-mer index 问题
- 2015年“深圳杯”数学建模夏令营-B题:DNA序列的k-mer index 问题
- 2015年“深圳杯”数学建模夏令营-B题:DNA序列的k-mer index 问题
- 数学建模有关DNA序列k-mer index的问题
- 数学建模_国2000A——DNA序列问题中的数据处理
- 数学建模常用Matlab/Lingo/c代码总结系列——最小费用最大流问题
- [刷题#1][FJOI2015]所有公共子序列问题
- 狱吏问题--数学建模的使用
- 数学建模(5)---煤矸石堆积问题
- 数学建模————统计问题之评价(三)
- 用动态规划法求解生物信息学中DNA序列比对的问题 (交叉学科应用实验)
- 数学建模-火箭发射问题
- Java解决数学建模之席位分配问题
- 数学建模常用Matlab/Lingo/c代码总结系列——整数规划问题
- 数学建模———线性规划问题
- 数学建模高级方法(十):时间序列模型
- 数学建模 TSP(旅行商问题) Lingo求解
- 数学建模(三)两个问题的思维过程
- 印刷厂问题(数学建模)
- 【数学建模】公平席位的分配问题