index.ts 5.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132
  1. import { spawn } from 'child_process';
  2. /* (c) Thomas Steimlé 2022
  3. * cat bwa_mem_splitters_on_HG38_Viral.sam | awk '$0~/^@/{next}{lxa=split($0,xa,"XA:Z:"); print $1"\t"$3"\t"$4; if(lxa>1){split(xa[2],xap,","); print $1"\t"xap[1]"\t"substr(xap[2],2)"\tXA"}}' | more
  4. * require os : cat, awk, sort, uniq
  5. *
  6. */
  7. const async_exec = (prog: string, args: string[], onData: Function, onErr: Function) => {
  8. return new Promise((resolve, reject) => {
  9. const child = spawn(prog, args, {shell: true})
  10. child.stdout.on('data', data => onData(data/*.toString().trim()*/))
  11. child.stderr.on('data', data => onErr(data.toString().trim()))
  12. child.on('error', err => reject(err))
  13. child.on('exit', code => resolve(code))
  14. })
  15. }
  16. const clusterSam = (
  17. input_sam: string | Array<string>,
  18. threshold: number,
  19. minReads : number
  20. ) => {
  21. return new Promise<any>( async (resolve, _reject) => {
  22. let inputSam: string = Array.isArray(input_sam) ? input_sam.join(' ') : input_sam
  23. let lineAcc: string = ''
  24. interface position {
  25. rname: string;
  26. position: number;
  27. }
  28. interface byContigs {
  29. [key: string]: position[]
  30. }
  31. let byContigs: byContigs = {}
  32. await async_exec('cat', [
  33. inputSam,
  34. '|',
  35. 'awk', '\'$0~/^@/{next}{lxa=split($0,xa,"XA:Z:"); print $1"\t"$3"\t"$4; if(lxa>1){split(xa[2],xap,","); print $1"\t"xap[1]"\t"substr(xap[2],2)"\tXA"}}\'', //skip header
  36. '|',
  37. 'sort',
  38. '|',
  39. 'uniq'
  40. ], (m: string) => {
  41. let tmpSeq: string[] = (lineAcc + m).split(/\n/)
  42. lineAcc = tmpSeq.pop() ! // 'uck typescript
  43. tmpSeq.map(e => {
  44. let tmpName: string = ''
  45. let tmpPos: position = {rname: '', position: 0}
  46. e.split(/\t/).map((el, i) => {
  47. switch (i) {
  48. case 0:
  49. tmpPos['rname'] = el
  50. break;
  51. case 1:
  52. tmpName = el
  53. break;
  54. case 2:
  55. tmpPos['position'] = Number(el)
  56. break;
  57. default:
  58. break;
  59. }
  60. })
  61. if (Array.isArray(byContigs[tmpName])) {
  62. byContigs[tmpName].push(tmpPos)
  63. } else {
  64. byContigs[tmpName] = [tmpPos]
  65. }
  66. })
  67. }, console.log)
  68. interface byReads {
  69. [key: string]: string[]
  70. }
  71. let byReads: byReads = {}
  72. interface posObj {
  73. [key: string]: string
  74. }
  75. interface posAll {
  76. [key: string]: posObj
  77. }
  78. let posAll: posAll = {}
  79. Object
  80. .keys(byContigs)
  81. .map(name => {
  82. let cluster = 0
  83. let firstPos = 0
  84. byContigs[name]
  85. .sort((a, b) => a.position - b.position)
  86. .map((e, i, a) => {
  87. if(i === 0) {
  88. firstPos = e.position
  89. }
  90. if (Math.abs(e.position - a[i-1]?.position) > threshold) {
  91. posAll[name][String(cluster)] = firstPos + '-' + a[i-1]?.position
  92. cluster = cluster + 1
  93. firstPos = e.position
  94. }
  95. // cluster = Math.abs(e.position - a[i-1]?.position) > threshold ? cluster + 1 : cluster
  96. const clutserName = cluster + '@' + name
  97. byReads[e.rname] = Array.isArray(byReads[e.rname]) ? [... new Set([...byReads[e.rname], clutserName])] : [clutserName]
  98. })
  99. })
  100. interface byClusters {
  101. [key: string]: string[]
  102. }
  103. let byClusters: byClusters = {}
  104. Object.keys(byReads).map(rname => {
  105. const tmpClusterName = byReads[rname].sort().map(e => {
  106. const splited = e.split(/@/)
  107. return splited[1] + ':' + posAll[splited[1]][splited[0]] + '(' + splited[0] + ')'
  108. }).join('<--->')
  109. byClusters[tmpClusterName] = Array.isArray(byClusters[tmpClusterName]) ? [... new Set([...byClusters[tmpClusterName], rname])] : [rname]
  110. })
  111. Object.keys(byClusters).map(e => byClusters[e].length < minReads ? delete byClusters[e] : null);
  112. resolve((Object.keys(byClusters).map(clusterName => ({clusterName, rnames: byClusters[clusterName]})).sort((a:any,b:any) => b.rnames.length - a.rnames.length)) )
  113. })
  114. }
  115. export { clusterSam }
  116. /*
  117. (async () => {
  118. console.log(await clusterSam('/home/thomas/Documents/Programmes/ttest/bwa_mem_splitters_on_HG38_Viral.sam', 333, 55));
  119. })()
  120. */