062add94-a5b5-4fd4-a562-e213ecc8da7c

Committed 20 Jan 2025 05:00PM UTC coverage: 86.049% (-0.1%) from 86.156%

Build # 062add94-a5b5-4fd4-a562-e213ecc8da7c

Build Type

push

circleci

Committed by

web-flow

Commit Message

Merge pull request #377 from DEploid-dev/shajoezhu-patch-1

Update config.yml

Run Details

5545 of 6444 relevant lines covered (86.05%)

5298085.9 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

50.43

/src/export/writeMcmcRelated.cpp

/*
 * dEploid is used for deconvoluting Plasmodium falciparum genome from
 * mix-infected patient sample.
 *
 * Copyright (C) 2016-2017 University of Oxford
 *
 * Author: Sha (Joe) Zhu
 *
 * This file is part of dEploid.
 *
 * dEploid is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 *
 */

#include "dEploidIO.hpp"
#include "mcmc.hpp"

void DEploidIO::writeMcmcRelated(McmcSample * mcmcSpl,
    string jobbrief, bool useIBD) {
    this->writeProp(mcmcSpl, jobbrief);
    this->writeLLK(mcmcSpl, jobbrief);
    this->writeHap(mcmcSpl->hap, jobbrief);

    if (useIBD == false) {
        this->writeVcf(mcmcSpl->hap,
            mcmcSpl->proportion.back(), jobbrief);
        this->siteOfTwoSwitchOne = mcmcSpl->siteOfTwoSwitchOne;
        this->siteOfTwoMissCopyOne = mcmcSpl->siteOfTwoMissCopyOne;
        this->siteOfTwoSwitchTwo = mcmcSpl->siteOfTwoSwitchTwo;
        this->siteOfTwoMissCopyTwo = mcmcSpl->siteOfTwoMissCopyTwo;
        this->siteOfOneSwitchOne = mcmcSpl->siteOfOneSwitchOne;
        this->siteOfOneMissCopyOne = mcmcSpl->siteOfOneMissCopyOne;

        this->finalSiteOfTwoSwitchOne = mcmcSpl->currentsiteOfTwoSwitchOne;
        this->finalSiteOfTwoMissCopyOne = mcmcSpl->currentsiteOfTwoMissCopyOne;
        this->finalSiteOfTwoSwitchTwo = mcmcSpl->currentsiteOfTwoSwitchTwo;
        this->finalSiteOfTwoMissCopyTwo = mcmcSpl->currentsiteOfTwoMissCopyTwo;
        this->finalSiteOfOneSwitchOne = mcmcSpl->currentsiteOfOneSwitchOne;
        this->finalSiteOfOneMissCopyOne = mcmcSpl->currentsiteOfOneMissCopyOne;

        // this->writeEventCount( );
    } else {
        // this->IBDpathChangeAt = mcmcSpl->IBDpathChangeAt;
        // this->finalIBDpathChangeAt = mcmcSpl->currentIBDpathChangeAt;
    }
}


void DEploidIO::writeProp(McmcSample * mcmcSample, string jobbrief) {
    string strExport = this->prefix_ + "." + jobbrief + ".prop";
    remove(strExport.c_str());
    ofstreamExportTmp.open(strExport.c_str(),
        ios::out | ios::app | ios::binary);
    for (size_t i = 0; i < mcmcSample->proportion.size(); i++) {
        for (size_t ii = 0; ii < mcmcSample->proportion[i].size(); ii++) {
            ofstreamExportTmp << setw(10) << mcmcSample->proportion[i][ii];
            ofstreamExportTmp << ((ii < (mcmcSample->proportion[i].size()-1)) ?
                "\t" : "\n");
        }
    }
    ofstreamExportTmp.close();
}


void DEploidIO::writeLLK(McmcSample * mcmcSample, string jobbrief) {
    string strExport = this->prefix_ + "." + jobbrief + ".llk";
    remove(strExport.c_str());
    ofstreamExportTmp.open(strExport.c_str(),
        ios::out | ios::app | ios::binary);
    for (size_t i = 0; i < mcmcSample->sumLLKs.size(); i++) {
        ofstreamExportTmp << mcmcSample->moves[i] << "\t"
                          << mcmcSample->sumLLKs[i] << endl;
    }
    ofstreamExportTmp.close();
}


void DEploidIO::writeHap(vector < vector <double> > &hap, string jobbrief) {
    string strExport = this->prefix_ + "." + jobbrief + ".hap";
    remove(strExport.c_str());
    ofstreamExportTmp.open(strExport.c_str(),
        ios::out | ios::app | ios::binary);
    // HEADER
    ofstreamExportTmp << "CHROM" << "\t" << "POS" << "\t";;
    for (size_t ii = 0; ii < kStrain_.getValue(); ii++) {
        ofstreamExportTmp << "h" << (ii+1);
        ofstreamExportTmp << ((ii < (kStrain_.getValue()-1)) ? "\t" : "\n");
    }

    size_t siteIndex = 0;
    for (size_t chromI = 0; chromI < chrom_.size(); chromI++) {
        dout << "chrom " << chromI << " length "
             << position_[chromI].size() << endl;
        for (size_t posI = 0; posI < position_[chromI].size(); posI++) {
            ofstreamExportTmp << chrom_[chromI]
                  << "\t" << static_cast<int>(position_[chromI][posI]) << "\t";
            for (size_t ii = 0; ii < hap[siteIndex].size(); ii++) {
                ofstreamExportTmp << hap[siteIndex][ii];
                ofstreamExportTmp << ((ii < (hap[siteIndex].size()-1)) ?
                    "\t" : "\n");
            }
            siteIndex++;
        }
    }

    assert(siteIndex == hap.size());
    ofstreamExportTmp.close();
}


void DEploidIO::writePanel(Panel *panel, size_t chromI, vector <string> hdr) {
    string strExport = this->prefix_ + ".panel." + to_string(chromI);
    remove(strExport.c_str());

    ofstreamExportTmp.open(strExport.c_str(),
        ios::out | ios::app | ios::binary);

    // HEADER
    ofstreamExportTmp << "CHROM" << "\t" << "POS" << "\t";;
    for (size_t ii = 0; ii < panel->truePanelSize_; ii++) {
        ofstreamExportTmp << hdr[ii];
        ofstreamExportTmp << ((ii < (panel->truePanelSize_-1)) ? "\t" : "\n");
    }

    size_t siteIndex = 0;
    for (size_t posI = 0; posI < position_[chromI].size(); posI++) {
        ofstreamExportTmp << chrom_[chromI] << "\t"
                          << static_cast<int>(position_[chromI][posI]) << "\t";
        for (size_t ii = 0; ii < panel->content_[siteIndex].size(); ii++) {
            ofstreamExportTmp << panel->content_[siteIndex][ii];
            ofstreamExportTmp << ((ii < (panel->content_[siteIndex].size()-1)) ?
                "\t" : "\n");
        }
        siteIndex++;
    }

    assert(siteIndex == panel->content_.size());
    ofstreamExportTmp.close();
}



void DEploidIO::writeVcf(vector < vector <double> > &hap,
    vector <double> &prop,
    string jobbrief) {
    if (!doExportVcf() ) return;

    this->strExportVcf = this->prefix_ + "." + jobbrief + ".vcf";
    if (compressVcf()) {
        strExportVcf += ".gz";
    }
    remove(strExportVcf.c_str());

    ogzstream ogstreamExport;
    ostream * writeTo;
    if (compressVcf()) {
        ogstreamExport.open(strExportVcf.c_str(), ios::out);
        writeTo = &ogstreamExport;
    } else {
        ofstreamExportTmp.open(strExportVcf.c_str(),
            ios::out | ios::app | ios::binary);
        writeTo = &ofstreamExportTmp;
    }

    // VCF HEADER
    if (this->useVcf()) {
        for (auto const& headerLine : this->vcfReaderPtr_->headerLines) {
            (*writeTo) << headerLine << endl;
        }
    } else {
        (*writeTo) << "##fileformat=VCFv4.2" << endl;
    }
    // DEploid call
    (*writeTo) << "##DEploid call: dEploid ";
    for (string s : argv_) {
        (*writeTo) << s << " ";
    }
    (*writeTo) << endl;

    // Include proportions
    for (size_t ii = 0; ii < prop.size(); ii++) {
        (*writeTo) << "##Proportion of strain "
                   << (this->useVcf() ? this->vcfReaderPtr_->sampleName_ : "h")
                   << "." << (ii+1)
                   << "=" << prop[ii] << endl;
    }

    // HEADER
    (*writeTo) << "#CHROM" << "\t"
               << "POS"    << "\t"
               << "ID"     << "\t"
               << "REF"    << "\t"
               << "ALT"    << "\t"
               << "QUAL"   << "\t"
               << "FILTER" << "\t"
               << "INFO"   << "\t"
               << "FORMAT" << "\t";
    for (size_t ii = 0; ii < kStrain_.getValue(); ii++) {
        (*writeTo) << (this->useVcf() ? this->vcfReaderPtr_->sampleName_ : "h")
                          << "." << (ii+1);
        (*writeTo) << ((ii < (kStrain_.getValue()-1)) ? "\t" : "\n");
    }

    size_t siteIndex = 0;
    for (size_t chromI = 0; chromI < chrom_.size(); chromI++) {
        for (size_t posI = 0; posI < position_[chromI].size(); posI++) {
            if (useVcf()) {
                (*writeTo) << this->vcfReaderPtr_->variants[siteIndex].chromStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].posStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].idStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].refStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].altStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].qualStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].filterStr
                   << "\t" << this->vcfReaderPtr_->variants[siteIndex].infoStr
                   << "\t" << "GT" << "\t";
            } else {
                (*writeTo) << chrom_[chromI]               << "\t"
                           << static_cast<int>(position_[chromI][posI]) << "\t"
                           << "."                          << "\t"
                           << "."                          << "\t"
                           << "."                          << "\t"
                           << "."                          << "\t"
                           << "."                          << "\t"
                           << "."                          << "\t"
                           << "GT"                         << "\t";
            }

            for (size_t ii = 0; ii < hap[siteIndex].size(); ii++) {
                (*writeTo) << hap[siteIndex][ii];
                (*writeTo) << ((ii < (hap[siteIndex].size()-1)) ? "\t" : "\n");
            }
            siteIndex++;
        }
    }

    assert(siteIndex == hap.size());
    if (compressVcf()) {
        ogstreamExport.close();
    } else {
        ofstreamExportTmp.close();
    }
}

1	/*
2	* dEploid is used for deconvoluting Plasmodium falciparum genome from
3	* mix-infected patient sample.
4	*
5	* Copyright (C) 2016-2017 University of Oxford
6	*
7	* Author: Sha (Joe) Zhu
8	*
9	* This file is part of dEploid.
10	*
11	* dEploid is free software: you can redistribute it and/or modify
12	* it under the terms of the GNU General Public License as published by
13	* the Free Software Foundation, either version 3 of the License, or
14	* (at your option) any later version.
15	*
16	* This program is distributed in the hope that it will be useful,
17	* but WITHOUT ANY WARRANTY; without even the implied warranty of
18	* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19	* GNU General Public License for more details.
20	*
21	* You should have received a copy of the GNU General Public License
22	* along with this program. If not, see <http://www.gnu.org/licenses/>.
23	*
24	*/
25
26	#include "dEploidIO.hpp"
27	#include "mcmc.hpp"
28
29	void DEploidIO::writeMcmcRelated(McmcSample * mcmcSpl,	4✔
30	string jobbrief, bool useIBD) {
31	this->writeProp(mcmcSpl, jobbrief);	8✔
32	this->writeLLK(mcmcSpl, jobbrief);	8✔
33	this->writeHap(mcmcSpl->hap, jobbrief);	4✔
34
35	if (useIBD == false) {	4✔
36	this->writeVcf(mcmcSpl->hap,	2✔
37	mcmcSpl->proportion.back(), jobbrief);
38	this->siteOfTwoSwitchOne = mcmcSpl->siteOfTwoSwitchOne;	2✔
39	this->siteOfTwoMissCopyOne = mcmcSpl->siteOfTwoMissCopyOne;	2✔
40	this->siteOfTwoSwitchTwo = mcmcSpl->siteOfTwoSwitchTwo;	2✔
41	this->siteOfTwoMissCopyTwo = mcmcSpl->siteOfTwoMissCopyTwo;	2✔
42	this->siteOfOneSwitchOne = mcmcSpl->siteOfOneSwitchOne;	2✔
43	this->siteOfOneMissCopyOne = mcmcSpl->siteOfOneMissCopyOne;	2✔
44
45	this->finalSiteOfTwoSwitchOne = mcmcSpl->currentsiteOfTwoSwitchOne;	2✔
46	this->finalSiteOfTwoMissCopyOne = mcmcSpl->currentsiteOfTwoMissCopyOne;	2✔
47	this->finalSiteOfTwoSwitchTwo = mcmcSpl->currentsiteOfTwoSwitchTwo;	2✔
48	this->finalSiteOfTwoMissCopyTwo = mcmcSpl->currentsiteOfTwoMissCopyTwo;	2✔
49	this->finalSiteOfOneSwitchOne = mcmcSpl->currentsiteOfOneSwitchOne;	2✔
50	this->finalSiteOfOneMissCopyOne = mcmcSpl->currentsiteOfOneMissCopyOne;	2✔
51
52	// this->writeEventCount( );
53	} else {
54	// this->IBDpathChangeAt = mcmcSpl->IBDpathChangeAt;
55	// this->finalIBDpathChangeAt = mcmcSpl->currentIBDpathChangeAt;
56	}
57	}	4✔
58
59
60	void DEploidIO::writeProp(McmcSample * mcmcSample, string jobbrief) {	4✔
61	string strExport = this->prefix_ + "." + jobbrief + ".prop";	8✔
62	remove(strExport.c_str());	4✔
63	ofstreamExportTmp.open(strExport.c_str(),	4✔
64	ios::out \| ios::app \| ios::binary);
65	for (size_t i = 0; i < mcmcSample->proportion.size(); i++) {	1,804✔
66	for (size_t ii = 0; ii < mcmcSample->proportion[i].size(); ii++) {	7,800✔
67	ofstreamExportTmp << setw(10) << mcmcSample->proportion[i][ii];	6,000✔
68	ofstreamExportTmp << ((ii < (mcmcSample->proportion[i].size()-1)) ?	6,000✔
69	"\t" : "\n");	7,800✔
70	}
71	}
72	ofstreamExportTmp.close();	4✔
73	}	4✔
74
75
76	void DEploidIO::writeLLK(McmcSample * mcmcSample, string jobbrief) {	4✔
77	string strExport = this->prefix_ + "." + jobbrief + ".llk";	8✔
78	remove(strExport.c_str());	4✔
79	ofstreamExportTmp.open(strExport.c_str(),	4✔
80	ios::out \| ios::app \| ios::binary);
81	for (size_t i = 0; i < mcmcSample->sumLLKs.size(); i++) {	1,804✔
82	ofstreamExportTmp << mcmcSample->moves[i] << "\t"	1,800✔
83	<< mcmcSample->sumLLKs[i] << endl;	1,800✔
84	}
85	ofstreamExportTmp.close();	4✔
86	}	4✔
87
88
89	void DEploidIO::writeHap(vector < vector <double> > &hap, string jobbrief) {	9✔
90	string strExport = this->prefix_ + "." + jobbrief + ".hap";	18✔
91	remove(strExport.c_str());	9✔
92	ofstreamExportTmp.open(strExport.c_str(),	9✔
93	ios::out \| ios::app \| ios::binary);
94	// HEADER
95	ofstreamExportTmp << "CHROM" << "\t" << "POS" << "\t";;	18✔
96	for (size_t ii = 0; ii < kStrain_.getValue(); ii++) {	38✔
97	ofstreamExportTmp << "h" << (ii+1);	29✔
98	ofstreamExportTmp << ((ii < (kStrain_.getValue()-1)) ? "\t" : "\n");	38✔
99	}
100
101	size_t siteIndex = 0;
102	for (size_t chromI = 0; chromI < chrom_.size(); chromI++) {	135✔
103	dout << "chrom " << chromI << " length "
104	<< position_[chromI].size() << endl;
105	for (size_t posI = 0; posI < position_[chromI].size(); posI++) {	4,735✔
106	ofstreamExportTmp << chrom_[chromI]
107	<< "\t" << static_cast<int>(position_[chromI][posI]) << "\t";	4,609✔
108	for (size_t ii = 0; ii < hap[siteIndex].size(); ii++) {	19,570✔
109	ofstreamExportTmp << hap[siteIndex][ii];	14,961✔
110	ofstreamExportTmp << ((ii < (hap[siteIndex].size()-1)) ?	14,961✔
111	"\t" : "\n");	19,570✔
112	}
113	siteIndex++;	4,609✔
114	}
115	}
116
117	assert(siteIndex == hap.size());
118	ofstreamExportTmp.close();	9✔
119	}	9✔
120
121
122	void DEploidIO::writePanel(Panel *panel, size_t chromI, vector <string> hdr) {	×
123	string strExport = this->prefix_ + ".panel." + to_string(chromI);	×
124	remove(strExport.c_str());	×
125
126	ofstreamExportTmp.open(strExport.c_str(),	×
127	ios::out \| ios::app \| ios::binary);
128
129	// HEADER
130	ofstreamExportTmp << "CHROM" << "\t" << "POS" << "\t";;	×
131	for (size_t ii = 0; ii < panel->truePanelSize_; ii++) {	×
132	ofstreamExportTmp << hdr[ii];
133	ofstreamExportTmp << ((ii < (panel->truePanelSize_-1)) ? "\t" : "\n");	×
134	}
135
136	size_t siteIndex = 0;
137	for (size_t posI = 0; posI < position_[chromI].size(); posI++) {	×
138	ofstreamExportTmp << chrom_[chromI] << "\t"
139	<< static_cast<int>(position_[chromI][posI]) << "\t";	×
140	for (size_t ii = 0; ii < panel->content_[siteIndex].size(); ii++) {	×
141	ofstreamExportTmp << panel->content_[siteIndex][ii];	×
142	ofstreamExportTmp << ((ii < (panel->content_[siteIndex].size()-1)) ?	×
143	"\t" : "\n");	×
144	}
145	siteIndex++;	×
146	}
147
148	assert(siteIndex == panel->content_.size());
149	ofstreamExportTmp.close();	×
150	}	×
151
152
153
154	void DEploidIO::writeVcf(vector < vector <double> > &hap,	3✔
155	vector <double> &prop,
156	string jobbrief) {
157	if (!doExportVcf() ) return;	3✔
158
159	this->strExportVcf = this->prefix_ + "." + jobbrief + ".vcf";	×
160	if (compressVcf()) {	×
161	strExportVcf += ".gz";
162	}
163	remove(strExportVcf.c_str());	×
164
165	ogzstream ogstreamExport;	×
166	ostream * writeTo;
167	if (compressVcf()) {	×
168	ogstreamExport.open(strExportVcf.c_str(), ios::out);
169	writeTo = &ogstreamExport;
170	} else {
171	ofstreamExportTmp.open(strExportVcf.c_str(),	×
172	ios::out \| ios::app \| ios::binary);
173	writeTo = &ofstreamExportTmp;	×
174	}
175
176	// VCF HEADER
177	if (this->useVcf()) {	×
178	for (auto const& headerLine : this->vcfReaderPtr_->headerLines) {	×
179	(*writeTo) << headerLine << endl;
180	}
181	} else {
182	(*writeTo) << "##fileformat=VCFv4.2" << endl;
183	}
184	// DEploid call
185	(*writeTo) << "##DEploid call: dEploid ";	×
186	for (string s : argv_) {	×
187	(*writeTo) << s << " ";	×
188	}
189	(*writeTo) << endl;
190
191	// Include proportions
192	for (size_t ii = 0; ii < prop.size(); ii++) {	×
193	(*writeTo) << "##Proportion of strain "
194	<< (this->useVcf() ? this->vcfReaderPtr_->sampleName_ : "h")	×
195	<< "." << (ii+1)	×
196	<< "=" << prop[ii] << endl;	×
197	}
198
199	// HEADER
200	(*writeTo) << "#CHROM" << "\t"
201	<< "POS" << "\t"
202	<< "ID" << "\t"
203	<< "REF" << "\t"
204	<< "ALT" << "\t"
205	<< "QUAL" << "\t"
206	<< "FILTER" << "\t"
207	<< "INFO" << "\t"
208	<< "FORMAT" << "\t";	×
209	for (size_t ii = 0; ii < kStrain_.getValue(); ii++) {	×
210	(*writeTo) << (this->useVcf() ? this->vcfReaderPtr_->sampleName_ : "h")	×
211	<< "." << (ii+1);	×
212	(*writeTo) << ((ii < (kStrain_.getValue()-1)) ? "\t" : "\n");	×
213	}
214
215	size_t siteIndex = 0;
216	for (size_t chromI = 0; chromI < chrom_.size(); chromI++) {	×
217	for (size_t posI = 0; posI < position_[chromI].size(); posI++) {	×
218	if (useVcf()) {	×
219	(*writeTo) << this->vcfReaderPtr_->variants[siteIndex].chromStr	×
220	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].posStr	×
221	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].idStr	×
222	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].refStr	×
223	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].altStr	×
224	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].qualStr	×
225	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].filterStr	×
226	<< "\t" << this->vcfReaderPtr_->variants[siteIndex].infoStr	×
227	<< "\t" << "GT" << "\t";	×
228	} else {
229	(*writeTo) << chrom_[chromI] << "\t"
230	<< static_cast<int>(position_[chromI][posI]) << "\t"
231	<< "." << "\t"
232	<< "." << "\t"
233	<< "." << "\t"
234	<< "." << "\t"
235	<< "." << "\t"
236	<< "." << "\t"
237	<< "GT" << "\t";	×
238	}
239
240	for (size_t ii = 0; ii < hap[siteIndex].size(); ii++) {	×
241	(*writeTo) << hap[siteIndex][ii];	×
242	(*writeTo) << ((ii < (hap[siteIndex].size()-1)) ? "\t" : "\n");	×
243	}
244	siteIndex++;	×
245	}
246	}
247
248	assert(siteIndex == hap.size());
249	if (compressVcf()) {	×
250	ogstreamExport.close();	×
251	} else {
252	ofstreamExportTmp.close();	×
253	}
254	}	×

DEploid-dev / DEploid / 062add94-a5b5-4fd4-a562-e213ecc8da7c

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous