<commit_before>#!/usr/bin/env node

var fs = require('fs'),
    path = require('path'),
    getopt = require('node-getopt'),
    MersenneTwister = require('mersennetwister'),
    assert = require('assert'),
    util = require('./util'),
    Ontology = require('./ontology'),
    Assocs = require('./assocs'),
    Model = require('./model'),
    MCMC = require('./mcmc'),
    Simulator = require('./simulator'),
    Benchmarker = require('./benchmarker')

var defaultSeed = 123456789
var defaultSamplesPerTerm = 100

// The default prior can be summarized as follows:
// - P(term present) = 1/#terms, sample size (#terms + 1)
// - P(false positive) = 1/#genes, sample size (#genes + 1)
// - P(false negative) = 1/#genes, sample size (#genes + 1)
var defaultTermPseudocount = 1
var defaultFalseNegPseudocount = 1
var defaultFalsePosPseudocount = 1

var defaultMoveRate = { flip: 1, swap: 1, randomize: 0 }
var defaultBenchReps = 1

var opt = getopt.create([
    ['o' , 'ontology=PATH'    , 'path to ontology file'],
    ['a' , 'assoc=PATH'       , 'path to gene-term association file'],
    ['g' , 'genes=PATH+'      , 'path to gene-set file(s)'],
    ['s' , 'samples=N'        , 'number of samples per term (default='+defaultSamplesPerTerm+')'],
    ['i' , 'ignore-missing'   , 'ignore missing terms & genes'],
    ['T',  'terms=N'          , 'pseudocount: active terms (default='+defaultTermPseudocount+')'],
    ['t',  'absent-terms=N'   , 'pseudocount: inactive terms (default=#terms)'],
    ['N',  'false-negatives=N', 'pseudocount: false negatives (default='+defaultFalseNegPseudocount+')'],
    ['p',  'true-positives=N' , 'pseudocount: true positives (default=#genes)'],
    ['P',  'false-positives=N', 'pseudocount: false positives (default='+defaultFalsePosPseudocount+')'],
    ['n',  'true-negatives=N' , 'pseudocount: true negatives (default=#genes)'],
    ['F',  'flip-rate=N'      , 'relative rate of term-toggling moves (default='+defaultMoveRate.flip+')'],
    ['S',  'swap-rate=N'      , 'relative rate of term-swapping moves (default='+defaultMoveRate.swap+')'],
    ['R',  'randomize-rate=N' , 'relative rate of term-randomizing moves (default='+defaultMoveRate.randomize+')'],
    ['l',  'log=TAG+'         , 'log various extra things (e.g. "move", "state")'],
    ['q' , 'quiet'            , 'don\'t log the usual things ("data", "progress")'],
    ['r' , 'rnd-seed=N'       , 'seed random number generator (default=' + defaultSeed + ')'],
    ['m' , 'simulate=N'       , 'instead of doing inference, simulate N gene sets'],
    ['x' , 'exclude-redundant', 'exclude redundant terms from simulation'],
    ['b' , 'benchmark'        , 'benchmark by running inference on simulated data'],
    ['B' , 'bench-reps=N'     , 'number of repetitions of benchmark (default='+defaultBenchReps+')'],
    ['h' , 'help'             , 'display this help message']
])              // create Getopt instance
.bindHelp()     // bind option 'help' to default action
.parseSystem(); // parse command line

function inputError(err) {
    throw new Error (err)
}

var ontologyPath = opt.options['ontology'] || inputError("You must specify an ontology")
var assocPath = opt.options['assoc'] || inputError("You must specify gene-term associations")

var samplesPerTerm = opt.options['samples'] ? parseInt(opt.options['samples']) : defaultSamplesPerTerm
var seed = opt.options['rnd-seed'] || defaultSeed

var logTags = opt.options['log'] || []
if (!opt.options['quiet'])
    logTags.push('data','progress')
function logging(tag) { return logTags.some(function(x) { return tag == x }) }

function readJsonFileSync (filename) {
    if (!fs.existsSync (filename))
        inputError ("File does not exist: " + filename)
    var data = fs.readFileSync (filename)
    return JSON.parse (data)
}

var ontologyJson = readJsonFileSync (ontologyPath)
var ontology = new Ontology ({termParents:ontologyJson})
if (logging('data'))
    console.warn("Read " + ontology.terms() + "-term ontology from " + ontologyPath)

var assocJson = readJsonFileSync (assocPath)
var assocs = new Assocs ({ ontology: ontology,
			   assocs: assocJson,
			   ignoreMissingTerms: opt.options['ignore-missing'] })

if (logging('data'))
    console.warn("Read " + assocs.nAssocs + " associations (" + assocs.genes() + " genes, " + assocs.relevantTerms().length + " terms) from " + assocPath)

var prior = {
    succ: {
	t: parseInt(opt.options['terms']) || defaultTermPseudocount,
	fp: parseInt(opt.options['false-positives']) || defaultFalsePosPseudocount,
	fn: parseInt(opt.options['false-negatives']) || defaultFalseNegPseudocount
    },
    fail: {
	t: parseInt(opt.options['absent-terms']) || assocs.relevantTerms().length,
	fp: parseInt(opt.options['true-negatives']) || assocs.genes(),
	fn: parseInt(opt.options['true-positives']) || assocs.genes()
    }
}

var moveRate = util.extend ({}, defaultMoveRate)
Object.keys(defaultMoveRate).forEach (function(r) {
    var arg = r + '-rate'
    if (arg in opt.options)
	moveRate[r] = parseInt (opt.options[arg])
})

var _generator
function generator() {
    _generator = _generator || new MersenneTwister (seed)
    return _generator
}

if (opt.options['benchmark'] || opt.options['bench-reps']) {
    var benchReps = opt.options['bench-reps'] ? parseInt(opt.options['bench-reps']) : defaultBenchReps
    var benchResults = { model: null, mcmc: null, benchmark: [] }
    var benchmarker = new Benchmarker ({ assocs: assocs })
    for (var benchRep = 0; benchRep < benchReps; ++benchRep) {

	if (logging('progress'))
	    console.warn("Starting benchmark repetition #" + (benchRep+1))

	var simResults = runSimulation()

	var genesJson = simResults.simulation.samples.map (function(sample) { return sample.gene.observed })
	if (logging('data'))
	    console.warn("Simulated " + genesJson.length + " gene set(s) of size [" + genesJson.map(function(l){return l.length}) + "]")

	var infResults = runInference (genesJson)

        benchmarker.add (simResults, infResults)
    }

    console.log(benchmarker)
    benchmarker.analyze()
    showResults (benchmarker.results)

} else if (opt.options['simulate']) {
    showResults (runSimulation())

} else {
    var genesPaths = opt.options['genes'] || inputError("You must specify at least one gene set")
    var genesJson = genesPaths.map (function(genesPath) { return readJsonFileSync (genesPath) })

    if (logging('data'))
	console.warn("Read " + genesJson.length + " gene set(s) of size [" + genesJson.map(function(l){return l.length}) + "] from [" + genesPaths + "]")

    var infResults = runInference (genesJson)
    showResults (infResults)
}

function runSimulation() {
    var sim = new Simulator ({ assocs: assocs,
			       generator: generator(),
			       prior: prior,
			       excludeRedundantTerms: opt.options['exclude-redundant'] })<commit_msg>refactored benchmark code ino separate file (removed debug message)<commit_after>#!/usr/bin/env node

var fs = require('fs'),
    path = require('path'),
    getopt = require('node-getopt'),
    MersenneTwister = require('mersennetwister'),
    assert = require('assert'),
    util = require('./util'),
    Ontology = require('./ontology'),
    Assocs = require('./assocs'),
    Model = require('./model'),
    MCMC = require('./mcmc'),
    Simulator = require('./simulator'),
    Benchmarker = require('./benchmarker')

var defaultSeed = 123456789
var defaultSamplesPerTerm = 100

// The default prior can be summarized as follows:
// - P(term present) = 1/#terms, sample size (#terms + 1)
// - P(false positive) = 1/#genes, sample size (#genes + 1)
// - P(false negative) = 1/#genes, sample size (#genes + 1)
var defaultTermPseudocount = 1
var defaultFalseNegPseudocount = 1
var defaultFalsePosPseudocount = 1

var defaultMoveRate = { flip: 1, swap: 1, randomize: 0 }
var defaultBenchReps = 1

var opt = getopt.create([
    ['o' , 'ontology=PATH'    , 'path to ontology file'],
    ['a' , 'assoc=PATH'       , 'path to gene-term association file'],
    ['g' , 'genes=PATH+'      , 'path to gene-set file(s)'],
    ['s' , 'samples=N'        , 'number of samples per term (default='+defaultSamplesPerTerm+')'],
    ['i' , 'ignore-missing'   , 'ignore missing terms & genes'],
    ['T',  'terms=N'          , 'pseudocount: active terms (default='+defaultTermPseudocount+')'],
    ['t',  'absent-terms=N'   , 'pseudocount: inactive terms (default=#terms)'],
    ['N',  'false-negatives=N', 'pseudocount: false negatives (default='+defaultFalseNegPseudocount+')'],
    ['p',  'true-positives=N' , 'pseudocount: true positives (default=#genes)'],
    ['P',  'false-positives=N', 'pseudocount: false positives (default='+defaultFalsePosPseudocount+')'],
    ['n',  'true-negatives=N' , 'pseudocount: true negatives (default=#genes)'],
    ['F',  'flip-rate=N'      , 'relative rate of term-toggling moves (default='+defaultMoveRate.flip+')'],
    ['S',  'swap-rate=N'      , 'relative rate of term-swapping moves (default='+defaultMoveRate.swap+')'],
    ['R',  'randomize-rate=N' , 'relative rate of term-randomizing moves (default='+defaultMoveRate.randomize+')'],
    ['l',  'log=TAG+'         , 'log various extra things (e.g. "move", "state")'],
    ['q' , 'quiet'            , 'don\'t log the usual things ("data", "progress")'],
    ['r' , 'rnd-seed=N'       , 'seed random number generator (default=' + defaultSeed + ')'],
    ['m' , 'simulate=N'       , 'instead of doing inference, simulate N gene sets'],
    ['x' , 'exclude-redundant', 'exclude redundant terms from simulation'],
    ['b' , 'benchmark'        , 'benchmark by running inference on simulated data'],
    ['B' , 'bench-reps=N'     , 'number of repetitions of benchmark (default='+defaultBenchReps+')'],
    ['h' , 'help'             , 'display this help message']
])              // create Getopt instance
.bindHelp()     // bind option 'help' to default action
.parseSystem(); // parse command line

function inputError(err) {
    throw new Error (err)
}

var ontologyPath = opt.options['ontology'] || inputError("You must specify an ontology")
var assocPath = opt.options['assoc'] || inputError("You must specify gene-term associations")

var samplesPerTerm = opt.options['samples'] ? parseInt(opt.options['samples']) : defaultSamplesPerTerm
var seed = opt.options['rnd-seed'] || defaultSeed

var logTags = opt.options['log'] || []
if (!opt.options['quiet'])
    logTags.push('data','progress')
function logging(tag) { return logTags.some(function(x) { return tag == x }) }

function readJsonFileSync (filename) {
    if (!fs.existsSync (filename))
        inputError ("File does not exist: " + filename)
    var data = fs.readFileSync (filename)
    return JSON.parse (data)
}

var ontologyJson = readJsonFileSync (ontologyPath)
var ontology = new Ontology ({termParents:ontologyJson})
if (logging('data'))
    console.warn("Read " + ontology.terms() + "-term ontology from " + ontologyPath)

var assocJson = readJsonFileSync (assocPath)
var assocs = new Assocs ({ ontology: ontology,
			   assocs: assocJson,
			   ignoreMissingTerms: opt.options['ignore-missing'] })

if (logging('data'))
    console.warn("Read " + assocs.nAssocs + " associations (" + assocs.genes() + " genes, " + assocs.relevantTerms().length + " terms) from " + assocPath)

var prior = {
    succ: {
	t: parseInt(opt.options['terms']) || defaultTermPseudocount,
	fp: parseInt(opt.options['false-positives']) || defaultFalsePosPseudocount,
	fn: parseInt(opt.options['false-negatives']) || defaultFalseNegPseudocount
    },
    fail: {
	t: parseInt(opt.options['absent-terms']) || assocs.relevantTerms().length,
	fp: parseInt(opt.options['true-negatives']) || assocs.genes(),
	fn: parseInt(opt.options['true-positives']) || assocs.genes()
    }
}

var moveRate = util.extend ({}, defaultMoveRate)
Object.keys(defaultMoveRate).forEach (function(r) {
    var arg = r + '-rate'
    if (arg in opt.options)
	moveRate[r] = parseInt (opt.options[arg])
})

var _generator
function generator() {
    _generator = _generator || new MersenneTwister (seed)
    return _generator
}

if (opt.options['benchmark'] || opt.options['bench-reps']) {
    var benchReps = opt.options['bench-reps'] ? parseInt(opt.options['bench-reps']) : defaultBenchReps
    var benchResults = { model: null, mcmc: null, benchmark: [] }
    var benchmarker = new Benchmarker ({ assocs: assocs })
    for (var benchRep = 0; benchRep < benchReps; ++benchRep) {

	if (logging('progress'))
	    console.warn("Starting benchmark repetition #" + (benchRep+1))

	var simResults = runSimulation()

	var genesJson = simResults.simulation.samples.map (function(sample) { return sample.gene.observed })
	if (logging('data'))
	    console.warn("Simulated " + genesJson.length + " gene set(s) of size [" + genesJson.map(function(l){return l.length}) + "]")

	var infResults = runInference (genesJson)

        benchmarker.add (simResults, infResults)
    }

    benchmarker.analyze()
    showResults (benchmarker.results)

} else if (opt.options['simulate']) {
    showResults (runSimulation())

} else {
    var genesPaths = opt.options['genes'] || inputError("You must specify at least one gene set")
    var genesJson = genesPaths.map (function(genesPath) { return readJsonFileSync (genesPath) })

    if (logging('data'))
	console.warn("Read " + genesJson.length + " gene set(s) of size [" + genesJson.map(function(l){return l.length}) + "] from [" + genesPaths + "]")

    var infResults = runInference (genesJson)
    showResults (infResults)
}

function runSimulation() {
    var sim = new Simulator ({ assocs: assocs,
			       generator: generator(),
			       prior: prior,
			       excludeRedundantTerms: opt.options['exclude-redundant'] })
{'fraction_non_alphanumeric': 0.16480140443274083, 'fraction_numerical': 0.003949967083607637, 'mean_word_length': 4.013568023469014, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2582064', 'n_tokens_mistral': 4536, 'n_tokens_neox': 4248, 'n_words': 1073}