地図情報作成・プログラム本体

方針を示します。

IMG_7236.HEIC

1039127 Bytes

ファイルダウンロードについて

ダウンロード

まず初めに、酵素１個を用いた時のアッセイの結果があるので、ここから絞れるその酵素の切断サイトの候補をvector<vector<ll>>に収納します。問題となる絞り方ですが、あるサイトにある・ないで場合分けすると到底間に合わないので、少し工夫をします。例えば上のように、挿入前の配列が30bp、挿入部位が20bpでアッセイにより2,3,4,41bpのバンドが出てきたとします。41は必ず挿入前の配列を全部含むので（問題の条件より、制限酵素の切断サイトは挿入前の配列ではなく挿入部位の配列にのみ存在するから）、2,3,4となるフラグメントの、挿入部位中の位置を考えれば良いことになります。上図のように塩基と切断サイトに番号を振っているので比較的考えやすいです。2,3,4のフラグメントはセットとしてあるので、まず合計9bpの塊が挿入部位中のどこにあるかで数え上げをします。次に、2,3,4の順列を考えてあげることで、考えられる候補を列挙して収納することができます。next_permutationをうまく使えばこの辺はできます。

ここで得られた候補を使い、その次の二つの酵素を使った結果と一致するかを調べます。先程の候補のやつを全部回せばこれは容易にできます。下のコードでは、酵素A,Bの候補の位置情報を取り出した後、vectorの機能を用いて重複をなくした後、そこから導出できる、フラグメントの長さ情報の理論値を入力にある実験値と比較し、一致しなければ他のを調べるということにしています。

では、コードを下に書いておきます(c++)。

#include <bits/stdc++.h>
using namespace std;
#define rep(i,n) for(ll i=0;i<n;i++)
#define repl(i,l,r) for(ll i=(l);i<(r);i++)
#define all(x) (x).begin(),(x).end()
using ll=long long;
using ld=long double;
using vl=vector<ll>;
using vvl=vector<vector<ll>>;
using pl=pair<ll,ll>;
using vpl=vector<pl>;
ll pow(ll a, ll b) {
ll x=1;
rep(i,b) {
x=x*a;
}
return x;
}
ll fact(ll n) {
ll a=1;
rep(i,n) {
a=a*(n-i);
}
return a;
}
int main() {
ll N;
ll M;
ll L;
cin>>N>>M>>L;
vvl Single(L);
vvl Double((L*(L-1))/2);
rep(i,L) {
ll k;
cin>>k;
rep(j,k) {
ll a;
cin>>a;
Single.at(i).push_back(a);
}
}

// candidates for cleaving sites of each enzyme are stored by using vector<vector<long long>>
vvl vec0(1000); //for enzyme0
vvl vec1(1000); //for enzyme 1
vvl vec2(1000); //for enzyme 2
vvl vec3(1000); // for enzyme 3

//Firstly Create lists of possible patterns of cleaving sites for each enzyme by using information gathered by experiments in which only one enzyme is used
ll count0=0;
ll count1=0;
ll count2=0;
ll count3=0;
ll subcount0=0;
ll subcount1=0;
ll subcount2=0;
ll subcount3=0;


rep(i,4) {
ll a1=Single.at(i).size()-1;
ll a[109];
ll b[fact(a1)][a1];
for (ll i = 0; i < a1; i++) {
a[i]=i+1;
}
ll count=0;

do {
for (ll i = 0; i < a1; i++) {
if (i+1)


b[count][i]=a[i]-1;
}

count++;
} while(next_permutation(a, a + a1));

if(i==1) {
ll count=0;
subcount1=a1+1;
rep(j,a1) {
count+=Single.at(i).at(j);
}
ll x=M-count+1;
ll c=0;
rep(j,x) {
rep(k,fact(a1)) {

ll keep=N-1+j;
vec1.at(c).push_back(keep);
count1++;
ll d=0;
rep(m,a1) {
d+=Single.at(1).at(b[k][m]);

ll z=d+keep;

vec1.at(c).push_back(z);
count1++;

}
vec1.at(c).push_back(keep+count);

c++;

}//repk fact end

}//repj end

}//if i==1 end

else if(i==2) {
ll count=0;
subcount2=a1+1;
rep(j,a1) {
count+=Single.at(i).at(j);
}
ll x=M-count+1;
ll c=0;
rep(j,x) {
rep(k,fact(a1)) {

ll keep=N-1+j;
vec2.at(c).push_back(keep);
count2++;
ll d=0;
rep(m,a1) {
d+=Single.at(2).at(b[k][m]);

ll z=d+keep;

vec2.at(c).push_back(z);
count2++;

}
vec2.at(c).push_back(keep+count);

c++;

}//repk fact end

}//repj end

}//if i==2 end

else if(i==3) {
ll count=0;
subcount3=a1+1;
rep(j,a1) {
count+=Single.at(i).at(j);
}
ll x=M-count+1;
ll c=0;
rep(j,x) {
rep(k,fact(a1)) {

ll keep=N-1+j;
vec3.at(c).push_back(keep);
count3++;
ll d=0;
rep(m,a1) {
d+=Single.at(3).at(b[k][m]);

ll z=d+keep;

vec3.at(c).push_back(z);
count3++;

}
vec3.at(c).push_back(keep+count);

c++;

}//repk fact end

}//repj end

}//if i==3 end

else if(i==0) {
ll count=0;
subcount0=a1+1;
rep(j,a1) {
count+=Single.at(i).at(j);
}
ll x=M-count+1;
ll c=0;
rep(j,x) {
rep(k,fact(a1)) {

ll keep=N-1+j;
vec0.at(c).push_back(keep);
count0++;
ll d=0;
rep(m,a1) {
d+=Single.at(0).at(b[k][m]);

ll z=d+keep;

vec0.at(c).push_back(z);
count0++;

}
vec0.at(c).push_back(keep+count);

c++;

}//repk fact end

}//repj end

}//if i==0 end

}//end rep0-4

count0=count0/subcount0;
count1=count1/subcount1;
count2=count2/subcount2;
count3=count3/subcount3;

//end of part 1
//Then by using the list, test if a particular pattern matches the following results, gained by experiments using four enzymes

vl v0;
vl v1;
vl v2;
vl v3;
vl v4;
vl v5;
rep(i,6) {
ll a;
cin>>a;
rep(j,a) {
ll b;
cin>>b;
if(i==0) {
v0.push_back(b);
}
else if(i==1) {
v1.push_back(b);
}
else if(i==2) {
v2.push_back(b);
}
else if(i==3) {
v3.push_back(b);
}
else if(i==4) {
v4.push_back(b);
}
else if(i==5) {
v5.push_back(b);
}

}

}

sort(v0.begin(),v0.end());
sort(v1.begin(),v1.end());
sort(v2.begin(),v2.end());
sort(v3.begin(),v3.end());
sort(v4.begin(),v4.end());
sort(v5.begin(),v5.end());

//confirmation
rep(i,count0) {
rep(j,count1) {
rep(k,count2) {
rep(l,count3) {


// site accomodation
vl ec0;
vl ec1;
vl ec2;
vl ec3;
rep(ii,subcount0) {
ec0.push_back(vec0.at(i).at(ii));
}
rep(jj,subcount1) {
ec1.push_back(vec1.at(j).at(jj));
}
rep(kk,subcount2) {
ec2.push_back(vec2.at(k).at(kk));
}
rep(ln,subcount3) {
ec3.push_back(vec3.at(l).at(ln));
}

//01
vl vec010=ec0;
vl vec011=ec1;
copy(vec011.begin(),vec011.end(),std::back_inserter(vec010));///vec010 contains the whole
sort(vec010.begin(),vec010.end());
vec010.erase(unique(vec010.begin(), vec010.end()), vec010.end());

ll x0=vec010.size();
vl piece01;
ll sum01=0;
rep(ii,x0-1) {
piece01.push_back(vec010.at(ii+1)-vec010.at(ii));
sum01+=vec010.at(ii+1)-vec010.at(ii);
}

piece01.push_back(M+N-sum01);
sort(all(piece01));
ll par01=1;
ll fin01=v0.size();
if(x0!=fin01) {
par01=0;
x0=min(x0,fin01);
}
rep(ii,x0){
if(v0.at(ii)!=piece01.at(ii)) {
par01=0;
//break;
}
}



//02
vl vec020=ec0;
vl vec021=ec2;
copy(vec021.begin(),vec021.end(),std::back_inserter(vec020));///vec020 contains the whole
sort(vec020.begin(),vec020.end());
vec020.erase(unique(vec020.begin(), vec020.end()), vec020.end());

ll x1=vec020.size();
vl piece02;
ll sum02=0;
rep(ii,x1-1) {
piece02.push_back(vec020.at(ii+1)-vec020.at(ii));
sum02+=vec020.at(ii+1)-vec020.at(ii);
}

piece02.push_back(M+N-sum02);
sort(all(piece02));
ll par02=1;
ll fin02=v1.size();
if(x1!=fin02) {
par02=0;
x1=min(x1,fin02);
}
rep(ii,x1){
if(v1.at(ii)!=piece02.at(ii)) {
par02=0;
//break;
}
}


//03
vl vec030=ec0;
vl vec031=ec3;
copy(vec031.begin(),vec031.end(),std::back_inserter(vec030));///vec030 contains the whole
sort(vec030.begin(),vec030.end());
vec030.erase(unique(vec030.begin(), vec030.end()), vec030.end());

ll x2=vec030.size();
vl piece03;
ll sum03=0;
rep(ii,x2-1) {
piece03.push_back(vec030.at(ii+1)-vec030.at(ii));
sum03+=vec030.at(ii+1)-vec030.at(ii);
}

piece03.push_back(M+N-sum03);
sort(all(piece03));
//cout<<piece03.size();
ll par03=1;
ll fin03=v2.size();
if(x2!=fin03) {
par03=0;
x2=min(x2,fin03);
}
rep(ii,x2){
if(v2.at(ii)!=piece03.at(ii)) {
par03=0;
//break;
}
}


//12
vl vec120=ec1;
vl vec121=ec2;
copy(vec121.begin(),vec121.end(),std::back_inserter(vec120));///vec120 contains the whole
sort(vec120.begin(),vec120.end());
vec120.erase(unique(vec120.begin(), vec120.end()), vec120.end());

ll x3=vec120.size();
vl piece12;
ll sum12=0;
rep(ii,x3-1) {
piece12.push_back(vec120.at(ii+1)-vec120.at(ii));
sum12+=vec120.at(ii+1)-vec120.at(ii);
}

piece12.push_back(M+N-sum12);
sort(all(piece12));
ll par12=1;
ll fin12=v3.size();
if(x3!=fin12) {
par12=0;
x3=min(x3,fin12);
}
rep(ii,x3){
if(v3.at(ii)!=piece12.at(ii)) {
par12=0;
//break;
}
}


//13
vl vec130=ec1;
vl vec131=ec3;
copy(vec131.begin(),vec131.end(),std::back_inserter(vec130));///vec130 contains the whole
sort(vec130.begin(),vec130.end());
vec130.erase(unique(vec130.begin(), vec130.end()), vec130.end());

ll x4=vec130.size();
vl piece13;
ll sum13=0;
rep(ii,x4-1) {
piece13.push_back(vec130.at(ii+1)-vec130.at(ii));
sum13+=vec130.at(ii+1)-vec130.at(ii);
}

piece13.push_back(M+N-sum13);
sort(all(piece13));
ll par13=1;
ll fin13=v4.size();
if(x4!=fin13) {
par13=0;
x4=min(x4,fin13);
}
rep(ii,x4){
if(v4.at(ii)!=piece13.at(ii)) {
par13=0;
//break;
}
}


//23
vl vec230=ec2;
vl vec231=ec3;
copy(vec231.begin(),vec231.end(),std::back_inserter(vec230));///vec230 contains the whole
sort(vec230.begin(),vec230.end());
vec230.erase(unique(vec230.begin(), vec230.end()), vec230.end());


ll x5=vec230.size();
vl piece23;
ll sum23=0;
rep(ii,x5-1) {
piece23.push_back(vec230.at(ii+1)-vec230.at(ii));
sum23+=vec230.at(ii+1)-vec230.at(ii);
}

piece23.push_back(M+N-sum23);
sort(all(piece23));
ll par23=1;
ll fin23=v5.size();
if(x5!=fin23) {
par23=0;
x5=min(x5,fin23);
}
rep(ii,x5){

if(v5.at(ii)!=piece23.at(ii)) {
par23=0;

//break;
}
}


if((par01==1)&&(par02==1)&&(par03==1)&&(par12==1)&&(par13==1)&&(par23==1)) {
cout<<0;
rep(ii,subcount0) {
cout<<" "<<vec0.at(i).at(ii);
}
cout<<endl;
cout<<1;
rep(ii,subcount1) {
cout<<" "<<vec1.at(j).at(ii);
}
cout<<endl;
cout<<2;
rep(ii,subcount2) {
cout<<" "<<vec2.at(k).at(ii);
}
cout<<endl;
cout<<3;
rep(ii,subcount3) {
cout<<" "<<vec3.at(l).at(ii);
}
cout<<endl;



//break;
}//endif

}
}
}
}
}//end all