Mercurial > repos > timpalpant > java_genomics_toolkit

comparison java-genomics-toolkit/src/edu/unc/genomics/visualization/MatrixAligner.java @ 0:1daf3026d231

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Upload alpha version
author timpalpant
date Mon, 13 Feb 2012 21:55:55 -0500
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:1daf3026d231
1 package edu.unc.genomics.visualization;
2
3 import java.io.BufferedWriter;
4 import java.io.IOException;
5 import java.nio.charset.Charset;
6 import java.nio.file.Files;
7 import java.nio.file.Path;
8 import java.util.Arrays;
9 import java.util.Iterator;
10 import java.util.List;
11
12 import org.apache.commons.lang3.ArrayUtils;
13 import org.apache.log4j.Logger;
14 import org.broad.igv.bbfile.WigItem;
15
16 import com.beust.jcommander.Parameter;
17
18 import edu.unc.genomics.BedEntry;
19 import edu.unc.genomics.CommandLineTool;
20 import edu.unc.genomics.ReadablePathValidator;
21 import edu.unc.genomics.io.BedFile;
22 import edu.unc.genomics.io.WigFile;
23 import edu.unc.genomics.io.WigFileException;
24
25 public class MatrixAligner extends CommandLineTool {
26
27 private static final Logger log = Logger.getLogger(MatrixAligner.class);
28
29 @Parameter(names = {"-i", "--input"}, description = "Input file (Wig)", required = true)
30 public WigFile inputFile;
31 @Parameter(names = {"-l", "--loci"}, description = "Loci file (Bed)", required = true, validateWith = ReadablePathValidator.class)
32 public Path lociFile;
33 @Parameter(names = {"-m", "--max"}, description = "Truncate width (base pairs)")
34 public Integer maxWidth;
35 @Parameter(names = {"-o", "--output"}, description = "Output file (matrix2png format)", required = true)
36 public Path outputFile;
37
38 private List<BedEntry> loci;
39
40 @Override
41 public void run() throws IOException {
42 log.debug("Loading alignment intervals");
43 try (BedFile bed = new BedFile(lociFile)) {
44 loci = bed.loadAll();
45 }
46
47 // Compute the matrix dimensions
48 int leftMax = Integer.MIN_VALUE;
49 int rightMax = Integer.MIN_VALUE;
50 for (BedEntry entry : loci) {
51 int left = Math.abs(entry.getValue().intValue()-entry.getStart());
52 int right = Math.abs(entry.getValue().intValue()-entry.getStop());
53 if (left > leftMax) {
54 leftMax = left;
55 }
56 if (right > rightMax) {
57 rightMax = right;
58 }
59 }
60
61 int m = loci.size();
62 int n = leftMax + rightMax + 1;
63 int alignmentPoint = leftMax;
64 log.info("Intervals aligned into: " + m+"x"+n + " matrix");
65 log.info("Alignment point: " + alignmentPoint);
66
67 int leftBound = 0;
68 int rightBound = n-1;
69 if (maxWidth != null && maxWidth < n) {
70 log.info("Truncated to: " + m+"x"+maxWidth);
71 int leftAlignDistance = alignmentPoint;
72 int rightAlignDistance = n - alignmentPoint - 1;
73 int halfMax = maxWidth / 2;
74
75 if (halfMax < leftAlignDistance && halfMax < rightAlignDistance) {
76 leftBound = alignmentPoint - halfMax;
77 rightBound = alignmentPoint + halfMax;
78 } else {
79 if (leftAlignDistance <= rightAlignDistance) {
80 rightBound = maxWidth;
81 } else {
82 leftBound = n - maxWidth;
83 }
84 }
85 }
86
87 log.debug("Initializing output file");
88 int count = 0, skipped = 0;
89 try (BufferedWriter writer = Files.newBufferedWriter(outputFile, Charset.defaultCharset())) {
90 writer.write("ID");
91 for (int i = leftBound-alignmentPoint; i <= rightBound-alignmentPoint; i++) {
92 writer.write("\t"+i);
93 }
94 writer.newLine();
95
96 log.debug("Iterating over all intervals");
97 String[] row = new String[n];
98 for (BedEntry entry : loci) {
99 Iterator<WigItem> result = null;
100 try {
101 result = inputFile.query(entry);
102 } catch (WigFileException e) {
103 skipped++;
104 continue;
105 }
106
107 float[] data = WigFile.flattenData(result, entry.getStart(), entry.getStop());
108 // Reverse if on the crick strand
109 if (entry.isCrick()) {
110 ArrayUtils.reverse(data);
111 }
112
113 // Position the data in the matrix
114 // Locus alignment point (entry value) should be positioned over the matrix alignment point
115 int n1 = alignmentPoint - Math.abs(entry.getValue().intValue()-entry.getStart());
116 int n2 = alignmentPoint + Math.abs(entry.getValue().intValue()-entry.getStop());
117 assert data.length == n2-n1+1;
118
119 Arrays.fill(row, "-");
120 for (int i = 0; i < data.length; i++) {
121 if (!Float.isNaN(data[i])) {
122 row[n1+i] = String.valueOf(data[i]);
123 }
124 }
125
126 // Write to output
127 String id = ((entry.getId() == null) ? entry.getId() : "Row "+(count++));
128 writer.write(id);
129 for (int i = leftBound; i <= rightBound; i++) {
130 writer.write("\t"+row[i]);
131 }
132 writer.newLine();
133 }
134 }
135
136 inputFile.close();
137 log.info(count + " intervals processed");
138 log.info(skipped + " intervals skipped");
139 }
140
141 public static void main(String[] args) {
142 new MatrixAligner().instanceMain(args);
143 }
144
145 }