diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/lmm.cpp | 33 |
1 files changed, 18 insertions, 15 deletions
diff --git a/src/lmm.cpp b/src/lmm.cpp index 9a1d318..324e310 100644 --- a/src/lmm.cpp +++ b/src/lmm.cpp @@ -1879,7 +1879,6 @@ void LMM::Analyze(std::function< SnpNameValues(size_t) >& fetch_snp, // Results for LMM. class SUMSTAT2 { public: - string marker; double beta; // REML estimator for beta. double se; // SE for beta. double lambda_remle; // REML estimator for lambda. @@ -2004,8 +2003,7 @@ void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp, time_opt += (clock() - time_start) / (double(CLOCKS_PER_SEC) * 60.0); // Store summary data. - SUMSTAT2 st = {marker, beta, se, lambda_remle, lambda_mle, - p_wald, p_lrt, p_score, logl_H1}; + SUMSTAT2 st = {beta, se, lambda_remle, lambda_mle, p_wald, p_lrt, p_score, logl_H1}; sumstat2.push_back(st); } }; @@ -2033,8 +2031,8 @@ void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp, auto success = get<0>(tup); if (!success) continue; - auto marker = get<1>(tup); - auto chr = get<2>(tup); + // auto marker = get<1>(tup); + // auto chr = get<2>(tup); auto mpos = get<3>(tup); auto gs = get<4>(tup); // cout << t << " SNP: " << snp << endl; @@ -2109,7 +2107,7 @@ void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp, auto sumstats = [&] (SUMSTAT2 st) { outfile << scientific << setprecision(6); - outfile << st.marker << "\t"; + // outfile << st.marker << "\t"; if (a_mode != M_LMM2) { outfile << st.beta << "\t"; outfile << st.se << "\t"; @@ -2215,7 +2213,8 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval, string marker; uint8_t chr; - size_t pos; + uint64_t pos; + // size_t pos; vector<double> gs; if (success) { @@ -2224,21 +2223,25 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval, const float* gsbuf = reinterpret_cast<const float*>(value.data()); // unpack chr+pos from key - if (key.size() != 5) { + if (key.size() != 10) { + cerr << "key.size=" << key.size() << endl; throw std::runtime_error("Invalid key size"); } + // "S>L>L>" const uint8_t* data = reinterpret_cast<const uint8_t*>(key.data()); - chr = static_cast<uint8_t>(data[0]); - // Extract big-endian uint32 manually - pos = (static_cast<uint32_t>(data[1]) << 24) | - (static_cast<uint32_t>(data[2]) << 16) | - (static_cast<uint32_t>(data[3]) << 8) | - (static_cast<uint32_t>(data[4])); + chr = static_cast<uint8_t>(data[1]); + // Extract big-endian uint32 manually abcd -> dcba + pos = (static_cast<uint32_t>(data[2]) << 24) | + (static_cast<uint32_t>(data[4]) << 8) | + (static_cast<uint32_t>(data[3]) << 16) | + (static_cast<uint32_t>(data[5])); + auto pos2 = __builtin_bswap32(pos); string_view value2; marker_mdb.get(rtxn,key,value2); marker = string(value2); - // cout << static_cast<int>(chr) << ":" << pos << " " << marker << endl ; + // 1 rs13476251 174792257 + cout << static_cast<int>(chr) << ":" << pos2 << "|" << pos << " " << marker << endl ; auto size = num_floats; gs.reserve(size); |