MCPcopy Create free account
hub / github.com/BoevaLab/FREEC / calculateLogLikelyHoodNormalMixtureForBAFs

Function calculateLogLikelyHoodNormalMixtureForBAFs

src/myFunc.cpp:2151–2382  ·  view source on GitHub ↗

Source from the content-addressed store, hash-verified

2149}
2150
2151double calculateLogLikelyHoodNormalMixtureForBAFs(vector <float> X,vector <float> mu,float middleComponentMinW, bool isMuFixed, bool CompleteGenomicsData) {
2152 double LogLikelyHood_i = -INFINITY ;
2153 int numberOfStates = mu.size();
2154 int N = X.size();
2155 int iterationCount = 0;
2156 int maxIterationCount = 10000;
2157
2158 float maxSigma = 0.07; //used only if isMuFixed==TRUE
2159 float minOmega = 0.6; //used only if isMuFixed==TRUE; for the sum of the two components in case of CN>2
2160 float minMinOmega = 0.15; //used only if isMuFixed==TRUE; for one component
2161
2162 //START TMP
2163// std::ofstream file ("/bioinfo/users/vboeva/Desktop/Neuroblastome/ANALYSES_CNG/Valentina_analysis/freec/rowBAFS_cn16.txt");
2164// for (int t = 0; t<N; t++) {
2165// file << X[t]<<"\n";
2166// }
2167// file.close();
2168 //END TMP
2169
2170 if (CompleteGenomicsData) { //correct mu:
2171 for (int i=0;i<numberOfStates;i++) {
2172 if (mu[i]!=HOMOZYG_MEAN)
2173 mu[i]*=0.82;
2174 }
2175 }
2176
2177 vector <float> Omega;
2178 //set Omega_0:
2179 if (numberOfStates==2) {
2180 Omega.push_back(0.9);
2181 Omega.push_back(0.1);
2182 } else if (numberOfStates==3) {
2183 Omega.push_back(0.18);
2184 Omega.push_back(0.02);
2185 Omega.push_back(0.8);
2186 } else if (numberOfStates==4) {
2187 Omega.push_back(0.18);
2188 Omega.push_back(0.02);
2189 Omega.push_back(0.4);
2190 Omega.push_back(0.4);
2191 } else {
2192 cerr << "Warning: in fitting by a mixture of normals, the maximum number of compontents is 4. Using "<<numberOfStates <<" components\n";
2193 return NA;
2194 }
2195
2196 vector <float> sigma;
2197 //set Sigma_0:
2198 if (numberOfStates==2) {
2199 sigma.push_back(0.2);
2200 sigma.push_back(0.2);
2201 } else if (numberOfStates==3) {
2202 sigma.push_back(0.2);
2203 sigma.push_back(0.2);
2204 sigma.push_back(maxSigma);
2205 } else if (numberOfStates==4) {
2206 sigma.push_back(0.2);
2207 sigma.push_back(0.2);
2208 sigma.push_back(maxSigma);

Callers 1

getBAFinfoFunction · 0.85

Calls 1

Tested by

no test coverage detected