diff options
-rw-r--r-- | gn/runner/gemma.scm | 43 | ||||
-rwxr-xr-x | scripts/precompute/precompute-hits.scm | 6 | ||||
-rwxr-xr-x | scripts/precompute/run-gemma.scm | 21 |
3 files changed, 42 insertions, 28 deletions
diff --git a/gn/runner/gemma.scm b/gn/runner/gemma.scm index 131284a..4c1c374 100644 --- a/gn/runner/gemma.scm +++ b/gn/runner/gemma.scm @@ -10,16 +10,37 @@ #:use-module (rnrs base) #:export ( + write-pheno-file run-gemma )) +(define (write-pheno-file fn traits) + (define bxd-inds (geno-inds-bxd "BXD.json")) + (assert (= 235 (length bxd-inds))) + (display bxd-inds) + (call-with-output-file fn + (lambda (port) + (for-each + (lambda (ind) + (begin + (let* [(value (assoc-ref traits ind)) + (outvalue (if value + value + "NA"))] + (format #t "~s ~s" ind outvalue) + (newline) + (display outvalue port) + (newline port)))) + bxd-inds) + (close port) + )) +) + (define (run-gemma population data-id probeset-id probesetfreeze-id name trait-name traits) "Run gemma-wrapper to compute GRM and GWA. On failure the run will stop(!)" - (define bxd-inds (geno-inds-bxd "BXD.json")) (define ids (string-append (int-to-string data-id) "," (int-to-string probeset-id) "," (int-to-string probesetfreeze-id))) - (assert (= 235 (length bxd-inds))) (if name (display (string-append "WE HAVE OUR " population " DATASET " name " and trait " trait-name " for precompute!\n"))) (display ids) @@ -31,24 +52,6 @@ (pheno-fn (string-append tmpdir "/pheno.txt")) (k-json-fn (string-append tmpdir "/K.json")) (gwa-json-fn (string-append tmpdir "/GWA.json"))] - (call-with-output-file pheno-fn - (lambda (port) - (for-each (lambda (ind) - (begin - (let* [(value (assoc-ref traits ind)) - (outvalue (if value - value - "NA"))] - (display outvalue) - (newline) - (display outvalue port) - (newline port)))) - bxd-inds) - (close port) - )) - - ;; set up with ./.guix-shell -- guile -L . -s ./scripts/precompute/precompute-hits.scm - ;; ---- to start GEMMA precompute inside container ;; env TMPDIR=. LD_LIBRARY_PATH=$GUIX_ENVIRONMENT/lib/ guile -L . -s ./scripts/precompute/precompute-hits.scm ;; --- First we compute K - control output goes to K.json diff --git a/scripts/precompute/precompute-hits.scm b/scripts/precompute/precompute-hits.scm index a18b33d..203dfbb 100755 --- a/scripts/precompute/precompute-hits.scm +++ b/scripts/precompute/precompute-hits.scm @@ -1,8 +1,8 @@ #! Precompute hits -This is the obsolete script to compute GEMMA directly from the DB. We have split out -this step into fetching traits and running gemma without DB. This allows scaling up -on HPC etc. +This is the *obsolete* script to compute GEMMA directly from the DB. We have split out +this step into separate steps for (1) fetching traits and (2) running gemma without DB. +This allows scaling up on HPC etc. Run from base dir with diff --git a/scripts/precompute/run-gemma.scm b/scripts/precompute/run-gemma.scm index b14927a..36edf8b 100755 --- a/scripts/precompute/run-gemma.scm +++ b/scripts/precompute/run-gemma.scm @@ -5,9 +5,9 @@ GEMMA on those. Run from base dir with -guix shell guile guile-dbi -- guile -L ../.. -s run-gemma.scm --help + ~/opt/guix-pull/bin/guix shell guile guile-dbi guile-json -- guile -L . -e main -s scripts/precompute/run-gemma.scm 115476.json -and with some extra paths +and with some extra paths (for gemma) . .guix-shell ruby --expose=/home/wrk/iwrk/opensource/code/genetics/gemma-wrapper/=/gemma-wrapper --expose=/home/wrk/iwrk/opensource/code/genetics/gemma/=/gemma -- env TMPDIR=tmp guile -L . -s ./scripts/precompute/run-gemma.scm @@ -33,10 +33,21 @@ and with some extra paths (help (single-char #\h) (value #f)))) (options (getopt-long args option-spec)) (help-wanted (option-ref options 'help #f))] + (display "RUNNING") (if help-wanted (format #t "list-traits-to-compute writes JSON traits files from the GN DB Usage: list-traits-to-compute [options...] -h, --help Display this help -") - (display "gemma-run") -))) +")) + (let [(trait-name "115475")] + (call-with-input-file "115475.json" + (lambda (port) + (let* [(json (json->scm port)) + (dataset (assoc-ref json "traits")) + (dataset-name (car (car dataset))) + (traits (assoc-ref dataset dataset-name)) + ] + (display dataset) + (write-pheno-file "pheno.txt" traits) + ; (run-gemma "BXD" data-id probeset-id probesetfreeze-id name trait-name traits))) + )))))) |