#include "Agent.h" #include "Processor.h" #include "Settings.h" #include #include #include #include #include #include #include #include #include #include #define sscanf_s scanf boost::random::random_device rng; boost::random::uniform_real_distribution<> uniform(0,1); void Processor::run(std::string fileName){ if(fileName.size() > 1){ this->initSettings(fileName); } else{ this->initSettings("C:/Deafness modelling/Deaf2.0 20.11.18/INIT.txt"); //Settings::OUTPUT_FILE = "D:/model/Deaf2.0 by Sergey 06.08/OUTPUT/output.deaf25.txt"; } int initPopSizeM = Settings::INIT_POP_SIZE_M; int initPopSizeW = Settings::INIT_POP_SIZE_W; ///float deafAlleleRatio = Settings::DEAF_ALLELE_RATIO; float birthRateH = Settings::BIRTH_RATE_H; // среднее число потомков float birthRateD = Settings::BIRTH_RATE_D; // среднее число потомков int maxGenerations = Settings::MAX_GENERATIONS; // Подготовка файла выдачи std::cout<<"Initializing output file..\t"; std::ofstream outputFile(Settings::OUTPUT_FILE); outputFile<<"Generation\tTotalPopulation\tMen\tWomen\tDeaf"; outputFile<<"\tDD marriages\tDD children"; outputFile<<"\tDH marriages\tDH children"; outputFile<<"\tHH marriages\tHH children"; outputFile<<"\tDA frequency\tDA homozygotes frequency"; // TODO Допишите ещё данные, которые вы хотите чтобы были outputFile< > matrixCurrent; //std::vector > matrixFuture; std::vector men; std::vector women; std::vector newMen; std::vector newWomen; std::default_random_engine& generator = Settings::RANDOM_GENERATOR; //generator.seed(std::time(0)); //RandomDevice::RandomDevice(unsigned long n) : rand_seed(n), engine(n){ } // Создаём мужчин std::cout<<"Creating men..("< distrUI(0.f,1.f); if(distrUI(generator) < Settings::SPONTANEOUS_DEAF) isDeaf = true; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); men.push_back(a); } // Создали мужчин - слышащих гомозигот // Слышащие гетерозиготы (в них могут быть спонтанные) for(int i = 0; i < heterozygores; i++){ bool isDeaf = false; int alleleFather = 1; // TODO: Чуть аккуратнее int alleleMother = 0; // Спонтанная глухота - только для негенетических std::uniform_real_distribution distrUI(0.f,1.f); if(distrUI(generator) < Settings::SPONTANEOUS_DEAF) isDeaf = true; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); men.push_back(a); } // Создали мужчин - слышащих гетерозигот // Глухие гомозиготы for(int i = 0; i < deafHomozygotes; i++){ bool isDeaf = true; int alleleFather = 1; int alleleMother = 1; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); std::uniform_real_distribution distrUI(0.f,1.f); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); men.push_back(a); }// Создали мужчин - неслышащих гомозигот std::random_shuffle(men.begin(), men.end()); // Создали всех мужчин std::cout<<"Done\n"; std::cout<<"Creating women..\t"; // Создаём женщин deafHomozygotes = (int) (initPopSizeW*Settings::DEAF_HOMOZYGOTES); hearHomozygotes = (int) (initPopSizeW*Settings::HEAR_HOMOZYGOTES); heterozygores = initPopSizeW - deafHomozygotes - hearHomozygotes; // Слышащие гомозиготы (в них могут быть спонтанные) for(int i = 0; i < hearHomozygotes; i++){ bool isDeaf = false; int alleleFather = 0; int alleleMother = 0; // Спонтанная глухота - только для негенетических std::uniform_real_distribution distrUI(0.f,1.f); if(distrUI(generator) < Settings::SPONTANEOUS_DEAF) isDeaf = true; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); women.push_back(a); } // Создали женщин - слышащих гомозигот // Слышащие гетерозиготы (в них могут быть спонтанные) for(int i = 0; i < heterozygores; i++){ bool isDeaf = false; int alleleFather = 1; // TODO: Чуть аккуратнее int alleleMother = 0; // Спонтанная глухота - только для негенетических std::uniform_real_distribution distrUI(0.f,1.f); if(distrUI(generator) < Settings::SPONTANEOUS_DEAF) isDeaf = true; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); women.push_back(a); } // Создали женщин - слышащих гетерозигот // Глухие гомозиготы for(int i = 0; i < deafHomozygotes; i++){ bool isDeaf = true; int alleleFather = 1; int alleleMother = 1; float mean = isDeaf ? Settings::SOCIAL_MEAN_D : Settings::SOCIAL_MEAN_H; float sigma= isDeaf ? Settings::SOCIAL_VAR_D : Settings::SOCIAL_VAR_H; std::normal_distribution distN(mean, sigma); float socialP = distN(generator); std::uniform_real_distribution distrUI(0.f,1.f); float signLangLevel = isDeaf ? Settings::SIGN_LANG_DEAF : Settings::SIGN_LANG_HEAR; bool signL = distrUI(generator) < signLangLevel ? true: false; Agent* a = new Agent(alleleFather, alleleMother, isDeaf, signL, socialP); women.push_back(a); }// Создали женщин - неслышащих гомозигот std::random_shuffle(women.begin(), women.end()); // Создали женщин std::cout<<"Done\n"; std::cout<<"\nMain cycle\n"; // Основной цикл for(int gen = 0; gen < maxGenerations; gen++){ std::cout<<"Generation\t"< vect(women.size(),0.0f); //matrixCurrent.assign(men.size(),vect); std::cout<<"\tDone\n"; // Эмулятор разреженной матрицы - пара <индекс, оценка> std::cout<<"\tInitializing settledPairs\t"; std::uniform_real_distribution distrUI(0.f,1.f); std::cout<<"."; std::vector > > settledPairs(men.size()); std::cout<<"."; for(int i = 0; i < settledPairs.size(); i++){ float f = distrUI(generator); int sign = f < 0.5f ? -1 : 1; int candidatePairs = Settings::CANDIDATE_PAIRS_MEAN + (int)(sign*f*Settings::CANDIDATE_PAIRS_VAR); // TODO: Несколько топологий краёв int left = i - candidatePairs/2; int right = i + candidatePairs/2; // HACK: 2018/11/23 Для глухих модель Егора if(Settings::DEAF_COMMUNITY_MODEL == 1){ left = 0; right = women.size() - 1; } if(left < 0){ right -= left; left = 0; } if(right >= women.size()){ right = women.size() - 1; left = right - candidatePairs; } std::vector > pairs(candidatePairs); for(unsigned int p = 0; p < pairs.size(); p++){ pairs.at(p).first = left; // Зараз и оценили друг-друга //matrixCurrent.at(i).at(left) pairs.at(p).second = esteem(men.at(i),women.at(left)); left++; } settledPairs.at(i) = pairs; } // (END) Эмулятор разреженной матрицы std::cout<<"Done\n"; // Заполняем матрицу взаимных оцениваний // По всем мужчинам /*for(unsigned int m = 0; m < men.size(); m++){ // По всем женщинам for(unsigned int w = 0; w < women.size(); w++){ if(matrixCurrent.at(m).at(w) >= 0.0f) matrixCurrent.at(m).at(w) = esteem(men.at(m),women.at(w)); } } */ // Заполнили матрицу взаимных оцениваний /* std::cout<<"\tGenerating offsprings\t"; // Формируем окончательные пары, генерируем потомков // По всем мужчинам std::cout << "MENS" << std::endl; for (auto m = 0; m < settledPairs.size(); ++m) { std::cout << "MEN: " << m << ": "; for (auto w = 0; w < settledPairs[m].size(); ++w) { std::cout << settledPairs[m][w].first << " : " << settledPairs[m][w].second << "; "; } std::cout << std::endl; } */ for(int m = 0; m < men.size(); m++){ // Ищем наиболее подходящую женщину int bestW = 0; float bestScore = 0.0f; /* БЫЛО ДО РАЗРЕЖЕННЫХ МАТРИЦ for(int w = 0; w < women.size(); w++){ if(matrixCurrent.at(m).at(w) > bestScore){ bestScore = matrixCurrent.at(m).at(w); bestW = w; } } // Нашли **/ // Вывод матрицы! // С РАЗРЕЖЕННЫМИ for(unsigned int p = 0; p < settledPairs.at(m).size(); p++){ unsigned int w = settledPairs.at(m).at(p).first; if(women.at(w)->alreadyMarried) continue; if(settledPairs.at(m).at(p).second/*matrixCurrent.at(m).at(w)*/ > bestScore){ bestScore = settledPairs.at(m).at(p).second; //matrixCurrent.at(m).at(w); bestW = w; } } // (END) С РАЗРЕЖЕННЫМ if(bestScore <= 0.0f) continue; // Этому мужчине никто не подошёл // Формируем пару, рожаем детей int* mariage = &marriageHH; int* childrenType = &childrenHH; float meanChildren = birthRateH; if((men.at(m)->deaf() && !women.at(bestW)->deaf()) || (!men.at(m)->deaf() && women.at(bestW)->deaf())){ mariage = &marriageDH; childrenType = &childrenDH; } if(men.at(m)->deaf() && women.at(bestW)->deaf()){ mariage = &marriageDD; childrenType = &childrenDD; } if(men.at(m)->deaf() || women.at(bestW)->deaf()){ meanChildren = birthRateD; } boost::math::beta_distribution<> distrB(Settings::BETA_A,Settings::BETA_B); float correction = (Settings::BETA_A+Settings::BETA_B)/Settings::BETA_A*meanChildren; // FOR WINDOWS float fChildren = distrB(generator)*correction; // float fChildren = quantile(distrB, uniform(rng)) * correction; int children = (int) std::floor(fChildren+0.5f); (*mariage)++; (*childrenType) += children; // TODO: Добавить знание жестового языка если хотя бы один сибс - глухой + знает жестовый язык std::vector newSibs; // Для учёта распространения жестового языка в одной семье std::vector deafSibs;// Для подсчёта вновь рождённых глухих for(int i = 0; i < children; i++){ Agent* newAgent = new Agent(men.at(m), m, women.at(bestW), bestW); newSibs.push_back(newAgent); if(newAgent->deaf()){ deafSibs.push_back(newAgent); } // Добавляем в общую популяцию std::uniform_int_distribution distrUI(0,1); if(distrUI(generator) == 0){ newMen.push_back(newAgent); } else{ newWomen.push_back(newAgent); } } // Трансмиссия жестового языка // Случай 1 - Хотя бы один родитель глухой и кто-то знает жестовый if(men.at(m)->deaf() || women.at(bestW)->deaf()){ if(men.at(m)->knowsSign() || women.at(bestW)->knowsSign()){ for(auto i = 0; i < newSibs.size(); i++){ newSibs.at(i)->teachSignLanguage(); } } } // Случай 2 - Хотя бы один сибс глухой else if(deafSibs.size() > 0){ std::uniform_real_distribution distrUF(0.f,1.f); // а) Он получает случайный шанс bool sibsKnowSL = false; // Глухие случайно узнают жестовый for(auto i = 0; i < deafSibs.size(); i++){ if(distrUF(generator) <= Settings::SIGN_LANG_DEAF){ deafSibs.at(i)->teachSignLanguage(); sibsKnowSL = true; } } // Слышащие случайно узнают жестовый for(auto i = 0; i < newSibs.size(); i++){ if(!newSibs.at(i)->deaf() && distrUF(generator) <= Settings::SIGN_LANG_HEAR){ newSibs.at(i)->teachSignLanguage(); sibsKnowSL = true; } } // б) Если родители или сибсы знают жестовый if(men.at(m)->knowsSign() || women.at(bestW)->knowsSign() || sibsKnowSL){ for(auto i = 0; i < newSibs.size(); i++){ newSibs.at(i)->teachSignLanguage(); } } } // Случай 3 - Слышащий просто случайно узнаёт жестовый else { std::uniform_real_distribution distrUF(0.f,1.f); for(auto i = 0; i < newSibs.size(); i++){ if(distrUF(generator) <= Settings::SIGN_LANG_HEAR){ newSibs.at(i)->teachSignLanguage(); } } } // (END) Трансмиссия жестового языка // Делаем эту женщину недоступной для всех остальных women.at(bestW)->alreadyMarried = true; //for(int i = m; i < men.size(); i++){ // matrixCurrent.at(i).at(bestW) = -100.0f; //} } // Конец генерации потомков std::cout<<"Done\n"; std::cout<<"\tSaving model state\t"; // Сохранить текущее состояние модели int totalP = men.size()+women.size(); outputFile<deaf()) totalDeaf++; } for(auto i = 0; i < women.size(); i++){ if(women.at(i)->deaf()) totalDeaf++; } outputFile<<"\t"<getDeafAllelesCount(); deafAllelesTotal += deafAlleles; //int increase = deafAllelesTotal == 2 ? 1 : 0; deafHomoz += deafAlleles/2; } for(unsigned int i = 0; i< women.size(); i++){ int deafAlleles = women.at(i)->getDeafAllelesCount(); deafAllelesTotal += deafAlleles; //int increase = deafAllelesTotal == 2 ? 1 : 0; deafHomoz += deafAlleles/2; } outputFile<<"\t"<getFatherID() == woman->getFatherID() || man->getMotherID() == woman->getMotherID()) && woman->getMotherID() != -1){ return -50.f; } if(man->deaf() == woman->deaf()){ if(man->deaf()){ ans += 2*Settings::WEIGHT_PHENO_D; float sum1 = man->knowsSign() ? Settings::WEIGHT_SIGN_D : 0; float sum2 = woman->knowsSign() ? Settings::WEIGHT_SIGN_D : 0; ans += sum1+sum2; } else{ ans += 2*Settings::WEIGHT_PHENO_H; } } else{ if(man->knowsSign() && woman->knowsSign()){ ans += Settings::WEIGHT_SIGN_D + Settings::WEIGHT_SIGN_H; } } if(man->getSocialP() + woman->getSocialP() > ans) return 0.0f; return ans; } void Processor::initSettings(std::string fileName){ std::vector script; std::ifstream file(fileName.c_str()); std::string str; while(getline(file,str)){ if(str.length() > 2) script.push_back(str); } file.close(); // process Script unsigned int curLine; for(curLine = 0; curLine < script.size(); curLine++){ if(script.at(curLine).at(0) == '/' || script.at(curLine).at(0) == '#') continue; std::vector tokensVector; // Search for tokens //boost::split(tokensVector, script.at(curLine), boost::is_any_of("=") ); boost::split(tokensVector, script.at(curLine), std::bind2nd(std::equal_to(), '=')); std::string instr = tokensVector.at(0); boost::to_lower(instr); boost::trim(instr); // Инструкции if(instr == "generations" || instr == "iterations" ){ int iter = Settings::MAX_GENERATIONS; sscanf_s(tokensVector.at(1).c_str(),"%d", &iter); Settings::MAX_GENERATIONS = iter; } if(instr == "init_pop_men"){ int pop = Settings::INIT_POP_SIZE_M; sscanf_s(tokensVector.at(1).c_str(),"%d", &pop); Settings::INIT_POP_SIZE_M = pop; } if(instr == "init_pop_women"){ int pop = Settings::INIT_POP_SIZE_W; sscanf_s(tokensVector.at(1).c_str(),"%d", &pop); Settings::INIT_POP_SIZE_W = pop; } if(instr == "candidate_pairs_mean"){ int cand = Settings::CANDIDATE_PAIRS_MEAN; sscanf_s(tokensVector.at(1).c_str(),"%d", &cand); Settings::CANDIDATE_PAIRS_MEAN = cand; } if(instr == "candidate_pairs_var"){ int cand = Settings::CANDIDATE_PAIRS_VAR; sscanf_s(tokensVector.at(1).c_str(),"%d", &cand); Settings::CANDIDATE_PAIRS_VAR = cand; } else if(instr == "spontaneous_deaf"){ float deaf = Settings::SPONTANEOUS_DEAF; sscanf_s(tokensVector.at(1).c_str(),"%f", &deaf); Settings::SPONTANEOUS_DEAF = deaf; } else if(instr == "deaf_allele_ratio"){ float deaf = Settings::DEAF_ALLELE_RATIO; sscanf_s(tokensVector.at(1).c_str(),"%f", &deaf); Settings::DEAF_ALLELE_RATIO = deaf; } else if(instr == "deaf_homozygotes"){ float deaf = Settings::DEAF_HOMOZYGOTES; sscanf_s(tokensVector.at(1).c_str(),"%f", &deaf); Settings::DEAF_HOMOZYGOTES = deaf; } else if(instr == "hear_homozygotes"){ float hear = Settings::HEAR_HOMOZYGOTES; sscanf_s(tokensVector.at(1).c_str(),"%f", &hear); Settings::HEAR_HOMOZYGOTES = hear; } else if(instr == "birth_rate_hearing" || instr == "birth_rate_h"){ float rate = Settings::BIRTH_RATE_H; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::BIRTH_RATE_H = rate; } else if(instr == "birth_rate_deaf" || instr == "birth_rate_d"){ float rate = Settings::BIRTH_RATE_D; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::BIRTH_RATE_D = rate; } else if(instr == "beta_a"){ float rate = Settings::BETA_A; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::BETA_A = rate; } else if(instr == "beta_b"){ float rate = Settings::BETA_B; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::BETA_B = rate; } else if(instr == "social_mean_h" || instr == "social_mean_hearing"){ float rate = Settings::SOCIAL_MEAN_H; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SOCIAL_MEAN_H = rate; } else if(instr == "social_mean_d" || instr == "social_mean_deaf"){ float rate = Settings::SOCIAL_MEAN_D; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SOCIAL_MEAN_D = rate; } else if(instr == "social_var_h" || instr == "social_var_hearing"){ float rate = Settings::SOCIAL_VAR_H; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SOCIAL_VAR_H = rate; } else if(instr == "social_var_d" || instr == "social_var_deaf"){ float rate = Settings::SOCIAL_VAR_D; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SOCIAL_VAR_D = rate; } else if(instr == "weight_phenotype_hearing" || instr == "weight_pheno_h"){ float rate = Settings::WEIGHT_PHENO_H; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::WEIGHT_PHENO_H = rate; } else if(instr == "weight_phenotype_deaf" || instr == "weight_pheno_d"){ float rate = Settings::WEIGHT_PHENO_D; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::WEIGHT_PHENO_D = rate; } else if(instr == "weight_sign_h" || instr == "weight_sign_hearing"){ float rate = Settings::WEIGHT_SIGN_H; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::WEIGHT_SIGN_H = rate; } else if(instr == "weight_sign_d" || instr == "weight_sign_deaf"){ float rate = Settings::WEIGHT_SIGN_D; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::WEIGHT_SIGN_D = rate; } else if(instr == "sign_lang_d" || instr == "sign_lang_deaf"){ float rate = Settings::SIGN_LANG_DEAF; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SIGN_LANG_DEAF = rate; } else if(instr == "sign_lang_h" || instr == "sign_lang_hear"){ float rate = Settings::SIGN_LANG_HEAR; sscanf_s(tokensVector.at(1).c_str(),"%f", &rate); Settings::SIGN_LANG_HEAR = rate; } else if(instr == "output" || instr == "output_file"){ std::string file = tokensVector.at(1); boost::trim(file); Settings::OUTPUT_FILE = file; } /*if(instr == "fullstat"){ int flag = 0; sscanf_s(tokensVector.at(1).c_str(),"%d", &flag); if(flag != 0) Settings::WRITE_FULL_STATISTICS_TO_FILE = true; }*/ /*else if(instr == "birth" || instr == "birthrate" || instr == "birth_rate" ){ double rate = Settings::BirthRate; sscanf_s(tokensVector.at(1).c_str(),"%lg", &rate); Settings::BirthRate = rate; }*/ } // (END) for(curLine = 0; curLine < script.size(); curLine++) }