Skip to content
Commits on Source (2)
Show whitespace changes
Inline Side-by-side
example/example.log
View file @ bae52791
+-----------------------------+
| |
| Eagle v2.4 |
| December 13, 2017 |
| Eagle v2.4.1 |
| November 18, 2018 |
| Po-Ru Loh |
| |
+-----------------------------+
Copyright (C) 2015-2017 Harvard University.
Copyright (C) 2015-2018 Harvard University.
Distributed under the GNU GPLv3+ open source license.
Command line options:
......@@ -56,25 +56,25 @@ Setting --histFactor=1.00
BEGINNING STEP 1
Time for step 1: 0.282519
Time for step 1 MN^2: 0.0201082
Time for step 1: 0.347481
Time for step 1 MN^2: 0.0270685
Making hard calls (time: 0.0164621)
Making hard calls (time: 0.0181408)
BEGINNING STEP 2
BATCH 1 OF 1
Building hash tables
.................................................................. (time: 0.0826769)
.................................................................. (time: 0.065696)
Phasing samples 1-379
Time for phasing batch: 0.335781
Time for phasing batch: 0.309231
Making hard calls (time: 0.0172811)
Making hard calls (time: 0.0192151)
Time for step 2: 0.435747
Time for step 2 MN^2: 0.0391194
Time for step 2: 0.394149
Time for step 2 MN^2: 0.0408464
BEGINNING STEP 3 (PBWT ITERS)
......@@ -87,108 +87,108 @@ BEGINNING PBWT ITER 1
BATCH 1 OF 10
Phasing samples 1-37
Time for phasing batch: 1.13401
Time for phasing batch: 0.977095
BATCH 2 OF 10
Phasing samples 38-75
Time for phasing batch: 1.04853
Time for phasing batch: 0.898498
BATCH 3 OF 10
Phasing samples 76-113
Time for phasing batch: 1.18956
Time for phasing batch: 0.950794
BATCH 4 OF 10
Phasing samples 114-151
Time for phasing batch: 1.1353
Time for phasing batch: 0.937661
BATCH 5 OF 10
Phasing samples 152-189
Time for phasing batch: 1.01419
Time for phasing batch: 0.909736
BATCH 6 OF 10
Phasing samples 190-227
Time for phasing batch: 1.05172
Time for phasing batch: 0.924443
BATCH 7 OF 10
Phasing samples 228-265
Time for phasing batch: 0.99322
Time for phasing batch: 0.9072
BATCH 8 OF 10
Phasing samples 266-303
Time for phasing batch: 0.99847
Time for phasing batch: 0.926957
BATCH 9 OF 10
Phasing samples 304-341
Time for phasing batch: 0.998471
Time for phasing batch: 0.888218
BATCH 10 OF 10
Phasing samples 342-379
Time for phasing batch: 1.05052
Time for phasing batch: 0.90299
Time for PBWT iter 1: 10.6142
Time for PBWT iter 1: 9.22363
BEGINNING PBWT ITER 2
BATCH 1 OF 10
Phasing samples 1-37
Time for phasing batch: 1.76314
Time for phasing batch: 1.49903
BATCH 2 OF 10
Phasing samples 38-75
Time for phasing batch: 1.79145
Time for phasing batch: 1.4698
BATCH 3 OF 10
Phasing samples 76-113
Time for phasing batch: 1.8253
Time for phasing batch: 1.51891
BATCH 4 OF 10
Phasing samples 114-151
Time for phasing batch: 1.54563
Time for phasing batch: 1.45062
BATCH 5 OF 10
Phasing samples 152-189
Time for phasing batch: 1.54299
Time for phasing batch: 1.42069
BATCH 6 OF 10
Phasing samples 190-227
Time for phasing batch: 1.59936
Time for phasing batch: 1.4851
BATCH 7 OF 10
Phasing samples 228-265
Time for phasing batch: 1.55847
Time for phasing batch: 1.44716
BATCH 8 OF 10
Phasing samples 266-303
Time for phasing batch: 1.51538
Time for phasing batch: 1.49382
BATCH 9 OF 10
Phasing samples 304-341
Time for phasing batch: 1.58919
Time for phasing batch: 1.43647
BATCH 10 OF 10
Phasing samples 342-379
Time for phasing batch: 1.64308
Time for phasing batch: 1.46241
Time for PBWT iter 2: 16.3741
Time for PBWT iter 2: 14.6841
Writing .haps.gz and .sample output
Time for writing output: 0.281509
Total elapsed time for analysis = 28.1508 sec
Time for writing output: 0.221324
Total elapsed time for analysis = 24.9776 sec
example/example_ref.log
View file @ bae52791
+-----------------------------+
| |
| Eagle v2.4 |
| December 13, 2017 |
| Eagle v2.4.1 |
| November 18, 2018 |
| Po-Ru Loh |
| |
+-----------------------------+
Copyright (C) 2015-2017 Harvard University.
Copyright (C) 2015-2018 Harvard University.
Distributed under the GNU GPLv3+ open source license.
Command line options:
......@@ -18,7 +18,7 @@ Command line options:
--outPrefix=target.phased
Setting number of threads to 1
Warning: The index file is older than the data file: ref.bcf.csi
[W::hts_idx_load2] The index file is older than the data file: ref.bcf.csi
Reference samples: Nref = 169
Target samples: Ntarget = 8
......@@ -43,7 +43,7 @@ Average # SNPs per cM: 42
Number of <=(64-SNP, 1cM) segments: 9
Average # SNPs per segment: 47
Time for reading input: 6.18623 sec
Time for reading input: 5.09263 sec
Fraction of heterozygous genotypes: 0.178114
Typical span of default 100-het history length: 13.44 cM
......@@ -59,15 +59,15 @@ PHASING ITER 1 OF 1
Phasing target samples
................................................................................
Time for phasing iter 1: 0.144591
Time for phasing iter 1: 0.132139
Writing vcf.gz output to target.phased.vcf.gz
[W::vcf_parse] INFO 'AC' is not defined in the header, assuming Type=String
[W::vcf_parse] INFO 'AN' is not defined in the header, assuming Type=String
[W::vcf_parse] INFO 'DP' is not defined in the header, assuming Type=String
[W::vcf_parse] INFO 'AFR_AF' is not defined in the header, assuming Type=String
[W::vcf_parse] INFO 'EX_TARGET' is not defined in the header, assuming Type=String
Time for writing output: 0.0239
Total elapsed time for analysis = 6.35535 sec
Time for writing output: 0.0356269
Total elapsed time for analysis = 5.26099 sec
Mean phase confidence of each target individual:
ID PHASE_CONFIDENCE
......
example/example_vcf.log
View file @ bae52791
+-----------------------------+
| |
| Eagle v2.4 |
| December 13, 2017 |
| Eagle v2.4.1 |
| November 18, 2018 |
| Po-Ru Loh |
| |
+-----------------------------+
Copyright (C) 2015-2017 Harvard University.
Copyright (C) 2015-2018 Harvard University.
Distributed under the GNU GPLv3+ open source license.
Command line options:
......@@ -39,25 +39,25 @@ Setting --histFactor=1.00
BEGINNING STEP 1
Time for step 1: 0.267283
Time for step 1 MN^2: 0.0164427
Time for step 1: 0.263286
Time for step 1 MN^2: 0.0227702
Making hard calls (time: 0.0152462)
Making hard calls (time: 0.018074)
BEGINNING STEP 2
BATCH 1 OF 1
Building hash tables
.................................................................. (time: 0.0786691)
.................................................................. (time: 0.0552511)
Phasing samples 1-379
Time for phasing batch: 0.331454
Time for phasing batch: 0.292423
Making hard calls (time: 0.0164559)
Making hard calls (time: 0.0192568)
Time for step 2: 0.426586
Time for step 2 MN^2: 0.0375612
Time for step 2: 0.366938
Time for step 2 MN^2: 0.0381578
BEGINNING STEP 3 (PBWT ITERS)
......@@ -70,108 +70,108 @@ BEGINNING PBWT ITER 1
BATCH 1 OF 10
Phasing samples 1-37
Time for phasing batch: 1.0778
Time for phasing batch: 1.02615
BATCH 2 OF 10
Phasing samples 38-75
Time for phasing batch: 1.0477
Time for phasing batch: 0.964455
BATCH 3 OF 10
Phasing samples 76-113
Time for phasing batch: 1.11441
Time for phasing batch: 0.988999
BATCH 4 OF 10
Phasing samples 114-151
Time for phasing batch: 1.06495
Time for phasing batch: 0.978426
BATCH 5 OF 10
Phasing samples 152-189
Time for phasing batch: 1.07034
Time for phasing batch: 0.935981
BATCH 6 OF 10
Phasing samples 190-227
Time for phasing batch: 1.13122
Time for phasing batch: 0.941438
BATCH 7 OF 10
Phasing samples 228-265
Time for phasing batch: 1.07535
Time for phasing batch: 0.935864
BATCH 8 OF 10
Phasing samples 266-303
Time for phasing batch: 1.16462
Time for phasing batch: 0.926938
BATCH 9 OF 10
Phasing samples 304-341
Time for phasing batch: 1.15337
Time for phasing batch: 0.927558
BATCH 10 OF 10
Phasing samples 342-379
Time for phasing batch: 1.17788
Time for phasing batch: 0.932059
Time for PBWT iter 1: 11.0777
Time for PBWT iter 1: 9.55792
BEGINNING PBWT ITER 2
BATCH 1 OF 10
Phasing samples 1-37
Time for phasing batch: 1.87655
Time for phasing batch: 1.56033
BATCH 2 OF 10
Phasing samples 38-75
Time for phasing batch: 1.6772
Time for phasing batch: 1.49036
BATCH 3 OF 10
Phasing samples 76-113
Time for phasing batch: 1.67862
Time for phasing batch: 1.51976
BATCH 4 OF 10
Phasing samples 114-151
Time for phasing batch: 1.57155
Time for phasing batch: 1.48635
BATCH 5 OF 10
Phasing samples 152-189
Time for phasing batch: 1.57551
Time for phasing batch: 1.43754
BATCH 6 OF 10
Phasing samples 190-227
Time for phasing batch: 1.64418
Time for phasing batch: 1.50166
BATCH 7 OF 10
Phasing samples 228-265
Time for phasing batch: 1.60561
Time for phasing batch: 1.48078
BATCH 8 OF 10
Phasing samples 266-303
Time for phasing batch: 1.63349
Time for phasing batch: 1.45079
BATCH 9 OF 10
Phasing samples 304-341
Time for phasing batch: 1.78741
Time for phasing batch: 1.47995
BATCH 10 OF 10
Phasing samples 342-379
Time for phasing batch: 1.80977
Time for phasing batch: 1.42872
Time for PBWT iter 2: 16.86
Time for PBWT iter 2: 14.8363
Writing vcf.gz output to phased.vcf.gz
Time for writing output: 0.204556
Total elapsed time for analysis = 34.7048 sec
Time for writing output: 0.176335
Total elapsed time for analysis = 30.3964 sec
src/DipTreePBWT.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/DipTreePBWT.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/Eagle.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......@@ -36,6 +36,7 @@
#include "omp.h"
#include <htslib/thread_pool.h>
#include <htslib/vcf.h>
#include "Types.hpp"
......@@ -3235,6 +3236,12 @@ namespace EAGLE {
cerr << "ERROR: Could not write to file " << outFile << endl;
exit(1);
}
// Get a threadpool for writing equal to
// the maximum number of threads that have been configured.
htsThreadPool p = {hts_tpool_init(omp_get_max_threads()),
0};
hts_set_thread_pool(out, &p);
hts_set_thread_pool(htsTmp, &p);
bcf_hdr_t *hdr = bcf_hdr_read(htsTmp);
bcf_hdr_append_eagle_version(hdr, argc, argv);
......@@ -3353,21 +3360,28 @@ namespace EAGLE {
bcf_hdr_destroy(hdr);
hts_close(out);
hts_close(htsTmp);
hts_tpool_destroy(p.pool);
remove(tmpFile.c_str());
}
// write phased output in non-ref mode
// differences from the above (ref-mode) are as follows:
// - does not take noImpMissing arg
// - checks chrom
// - does not increment m64j outside the bpStart-bpEnd region
// - does not delete tmpFile (now vcfFile with original input)
void Eagle::writeVcfNonRef(const string &vcfFile, const string &outFile, int inputChrom,
int chromX, double bpStart, double bpEnd, const string &writeMode,
bool keepMissingPloidyX, int argc, char **argv) const {
bool noImpMissing, bool keepMissingPloidyX, int argc, char **argv)
const {
htsFile *htsIn = hts_open(vcfFile.c_str(), "r");
htsFile *out = hts_open(outFile.c_str(), writeMode.c_str());
// Get a threadpool for writing equal to
// the maximum number of threads that have been configured.
htsThreadPool p = {hts_tpool_init(omp_get_max_threads()),
0};
hts_set_thread_pool(out, &p);
hts_set_thread_pool(htsIn, &p);
bcf_hdr_t *hdr = bcf_hdr_read(htsIn);
bcf_hdr_append_eagle_version(hdr, argc, argv);
......@@ -3423,7 +3437,7 @@ namespace EAGLE {
else if (!missing && minIdx > 0) { // ALT1/ALT2 het => don't phase (shouldn't happen)
ptr[0] = ptr[1] = bcf_gt_missing;
}
else { // REF/ALT* het => phase as called by Eagle
else if (!missing || !noImpMissing) { // miss||REF/ALT* het => phase as called by Eagle
for (int j = 0; j < ploidy; j++) {
uint64 nTargetHap = 2*i + j;
int altIdx = missing ? 1 : maxIdx;
......@@ -3436,7 +3450,7 @@ namespace EAGLE {
}
else { // haploid
mostRecentPloidy[i] = 1;
if ( bcf_gt_is_missing(ptr[0]) ) { // missing allele
if ( bcf_gt_is_missing(ptr[0]) && !noImpMissing ) { // missing allele
int j = 0;
uint64 nTargetHap = 2*i + j;
int altIdx = 1;
......@@ -3465,6 +3479,7 @@ namespace EAGLE {
bcf_hdr_destroy(hdr);
hts_close(out);
hts_close(htsIn);
hts_tpool_destroy(p.pool);
}
void Eagle::makeHardCalls(uint64 n0start, uint64 n0end, uint seed) {
......
src/Eagle.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......@@ -191,7 +191,7 @@ namespace EAGLE {
int argc, char**argv) const;
void writeVcfNonRef(const std::string &vcfFile, const std::string &outFile, int inputChrom,
int chromX, double bpStart, double bpEnd, const std::string &writeMode,
bool keepMissingPloidyX, int argc, char **argv) const;
bool noImpMissing, bool keepMissingPloidyX, int argc, char **argv) const;
void makeHardCalls(uint64 n0start, uint64 n0end, uint seed);
void initRefIter(int refIter);
void buildHashTables(int iter, int batch, int seed);
......
src/EagleImpMiss.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/EagleMain.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......@@ -83,8 +83,6 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
uint64 Nref = vcfData.getNref(), Ntarget = vcfData.getNtarget();
int iters = params.pbwtIters;
if (params.noImpMissing)
iters = 1;
if (iters == 0) {
if (Ntarget < Nref/2)
iters = 1;
......@@ -95,6 +93,7 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
cout << endl << "Auto-selecting number of phasing iterations: setting --pbwtIters to "
<< iters << endl << endl;
}
bool internalImpMissing = !params.noImpMissing || iters > 1; // internally impute missing
cout << endl << "BEGINNING PHASING" << endl;
......@@ -110,7 +109,7 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
HapBitsT *hapBitsTptr = NULL;
HapHedge *hapHedgePtr = NULL;
#ifdef OLD_IMP_MISSING
if (!params.noImpMissing) {
if (internalImpMissing) {
cout << "Making HapHedge" << endl;
hapBitsTptr = new HapBitsT(eagle.getHaploBitsT(), 2*eagle.getNlib(2+iter),
eagle.getMseg64(), eagle.getMaskSnps64j());
......@@ -138,10 +137,10 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
}
confs[i-Nref] = eagle.runPBWT
(i, nF1, nF2, params.Kpbwt/(iter<iters?2:1), params.expectIBDcM, params.histFactor,
iter==iters, iter>1, !params.noImpMissing, params.usePS, conPSall[i-Nref],
iter==iters, iter>1, internalImpMissing, params.usePS, conPSall[i-Nref],
params.chrom==params.chromX);
#ifdef OLD_IMP_MISSING
if (!params.noImpMissing) {
if (internalImpMissing) {
Timer tim;
eagle.imputeMissing(*hapHedgePtr, i);
timeImpMissing += tim.update_time();
......@@ -155,7 +154,7 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
}
}
#ifdef OLD_IMP_MISSING
if (!params.noImpMissing) {
if (internalImpMissing) {
delete hapHedgePtr;
delete hapBitsTptr;
}
......@@ -177,7 +176,7 @@ void phaseWithRef(EagleParams &params, Timer &timer, double t0, int argc, char *
cout << "Time for phasing iter " << iter << " MN^2: " << timeMN2 / params.numThreads
<< endl;
#ifdef OLD_IMP_MISSING
else if (!params.noImpMissing)
else if (internalImpMissing)
cout << "Time for phasing iter " << iter << " impMissing: "
<< timeImpMissing / params.numThreads << endl;
#endif
......@@ -216,7 +215,7 @@ int main(int argc, char *argv[]) {
cout << " +-----------------------------+" << endl;
cout << endl;
cout << "Copyright (C) 2015-2017 Harvard University." << endl;
cout << "Copyright (C) 2015-2018 Harvard University." << endl;
cout << "Distributed under the GNU GPLv3+ open source license." << endl << endl;
//cout << "Boost version: " << BOOST_LIB_VERSION << endl;
......@@ -503,8 +502,8 @@ int main(int argc, char *argv[]) {
string outFile = params.outPrefix + "." + params.vcfOutSuffix;
cout << "Writing " << params.vcfOutSuffix << " output to " << outFile << endl;
eagle.writeVcfNonRef(params.vcfFile, outFile, params.chrom, params.chromX, params.bpStart,
params.bpEnd, params.vcfWriteMode, params.keepMissingPloidyX, argc,
argv);
params.bpEnd, params.vcfWriteMode, params.noImpMissing,
params.keepMissingPloidyX, argc, argv);
}
else {
cout << "Writing .haps.gz and .sample output" << endl; timer.update_time();
......
src/EaglePBWT.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/EagleParams.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......@@ -237,10 +237,6 @@ namespace EAGLE {
if (vm.count("vcfRef") || vm.count("vcfTarget") || vm.count("vcf")) { // VCF mode
if (vm.count("vcf")) { // non-ref mode
if (noImpMissing) {
cerr << "ERROR: --noImpMissing is only supported in ref-mode" << endl;
return false;
}
if (bpFlanking != 0) {
cerr << "ERROR: --bpFlanking is only supported in ref-mode" << endl;
return false;
......@@ -255,13 +251,14 @@ namespace EAGLE {
cerr << "ERROR: --vcfTarget must be specified in reference-based phasing mode" << endl;
return false;
}
if (pbwtIters > 1 && noImpMissing) {
cerr << "ERROR: --pbwtIters cannot be greater than 1 if --noImpMissing is set" << endl;
return false;
}
if (vcfRef.substr(vcfRef.length()-4) != string(".bcf")) {
cerr << "WARNING: --vcfRef does not end in '.bcf'; BCF input is fastest" << endl;
}
if (chromStr=="0" && (bpStart>0 || bpEnd<1e9)) {
cerr << "ERROR: --chrom is required when specifying --bpStart or --bpEnd in ref mode"
<< endl;
return false;
}
}
// vcf input checks for both ref and non-ref mode
......@@ -291,7 +288,7 @@ namespace EAGLE {
return false;
}
}
else { // non-ref mode
else { // non-ref mode; plink input
if (famFile.empty()) {
cerr << "ERROR: fam file must be specified either using --fam or --bfile"
<< endl;
......@@ -308,7 +305,7 @@ namespace EAGLE {
return false;
}
if (noImpMissing) {
cerr << "ERROR: --noImpMissing is only supported in ref-mode" << endl;
cerr << "ERROR: --noImpMissing is only supported with vcf input" << endl;
return false;
}
if (bpFlanking != 0) {
......
src/EagleParams.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/FileUtils.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/FileUtils.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/GenoData.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......@@ -24,8 +24,11 @@
#include <cstring>
#include <cmath>
#include <htslib/thread_pool.h>
#include <htslib/vcf.h>
#include "omp.h"
#include "Types.hpp"
#include "FileUtils.hpp"
#include "MemoryUtils.hpp"
......@@ -156,7 +159,7 @@ namespace EAGLE {
}
if (!ret.empty() &&
(snp.chrom < ret.back().chrom ||
(snp.chrom == ret.back().chrom && (snp.physpos <= ret.back().physpos ||
(snp.chrom == ret.back().chrom && (snp.physpos < ret.back().physpos ||
snp.genpos < ret.back().genpos)))) {
if (numOutOfOrder < 5) {
cerr << "WARNING: Out-of-order snp in bim file: " << bimFile << endl;
......@@ -645,6 +648,14 @@ namespace EAGLE {
cout << endl;
cout << "Total post-QC indivs: N = " << N << endl;
cout << "Total post-QC SNPs: M = " << M << endl;
if (N < 2) {
cerr << "ERROR: N>=2 individuals are required" << endl;
exit(1);
}
if (M == 0) {
cerr << "ERROR: At least 1 SNP is required" << endl;
exit(1);
}
cout << "MAF spectrum: " << endl;
const double mafBounds6[7] = {0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.500001};
......@@ -680,10 +691,15 @@ namespace EAGLE {
bool noMapCheck, double cMmax) {
htsFile *fin = hts_open(vcfFile.c_str(), "r");
if (fin == NULL) {
cerr << "ERROR: Could not open " << vcfFile << " for reading" << endl;
exit(1);
}
htsThreadPool p = {hts_tpool_init(omp_get_max_threads()),
0};
hts_set_thread_pool(fin, &p);
bcf_hdr_t *hdr = bcf_hdr_read(fin);
bcf1_t *rec = bcf_init1();
int mgt = 0, *gt = NULL;
......@@ -768,6 +784,7 @@ namespace EAGLE {
bcf_destroy(rec);
bcf_hdr_destroy(hdr);
hts_close(fin);
hts_tpool_destroy(p.pool);
cout << "Read M = " << snpsPreQC.size() << " variants" << endl;
processMap(snpsPreQC, geneticMapFile, noMapCheck); // modify snpsPreQC
......
src/GenoData.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/HapHedge.cpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......
src/HapHedge.hpp
View file @ bae52791
/*
This file is part of the Eagle haplotype phasing software package
developed by Po-Ru Loh. Copyright (C) 2015-2016 Harvard University.
developed by Po-Ru Loh. Copyright (C) 2015-2018 Harvard University.
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
......