001/*
002 *                    BioJava development code
003 *
004 * This code may be freely distributed and modified under the
005 * terms of the GNU Lesser General Public Licence.  This should
006 * be distributed with the code.  If you do not have a copy,
007 * see:
008 *
009 *      http://www.gnu.org/copyleft/lesser.html
010 *
011 * Copyright for this code is held jointly by the individual
012 * authors.  These should be listed in @author doc comments.
013 *
014 * For more information on the BioJava project and its aims,
015 * or to join the biojava-l mailing list, visit the home page
016 * at:
017 *
018 *      http://www.biojava.org/
019 *
020 */
021package org.biojava.nbio.survival.cox.matrix;
022
023
024/**
025 *
026 * http://introcs.cs.princeton.edu/java/22library/Matrix.java.html
027 */
028/**
029 * ***********************************************************************
030 * Compilation: javac Matrix.java Execution: java Matrix
031 *
032 * A bare-bones collection of static methods for manipulating matrices.
033 *
034 ************************************************************************
035 */
036public class Matrix {
037
038
039
040        /**
041         *
042         * @param A
043         * @return
044         */
045        public static double[][] sqrt(double[][] A) {
046                double[][] d = new double[A.length][A[0].length];
047                for (int i = 0; i < d.length; i++) {
048                        for (int j = 0; j < d[0].length; j++) {
049                                d[i][j] = Math.sqrt(A[i][j]);
050                        }
051                }
052                return d;
053        }
054
055        /**
056         *
057         * @param A
058         * @return
059         */
060        public static double[] sqrt(double[] A) {
061                double[] d = new double[A.length];
062                for (int i = 0; i < d.length; i++) {
063                        d[i] = Math.sqrt(A[i]);
064                }
065                return d;
066        }
067
068        /**
069         *
070         * @param A
071         * @return
072         */
073        public static double[][] oneDivide(double[][] A) {
074                double[][] d = new double[A.length][A[0].length];
075                for (int i = 0; i < d.length; i++) {
076                        for (int j = 0; j < d[0].length; j++) {
077                                d[i][j] = 1.0 / A[i][j];
078                        }
079                }
080                return d;
081        }
082
083        /**
084         *
085         * @param A
086         * @return
087         */
088        public static double[] oneDivide(double[] A) {
089                double[] d = new double[A.length];
090                for (int i = 0; i < d.length; i++) {
091                        d[i] = 1.0 / A[i];
092                }
093                return d;
094        }
095
096        /**
097         *
098         * @param A
099         * @return
100         */
101        public static double[] diag(double[][] A) {
102                double[] d = new double[A.length];
103                for (int i = 0; i < d.length; i++) {
104                        d[i] = A[i][i];
105                }
106                return d;
107        }
108
109        /**
110         *
111         * @param A
112         * @return
113         */
114        public static double[][] diag(double[] A) {
115                double[][] d = new double[A.length][A.length];
116                for (int i = 0; i < d.length; i++) {
117                        d[i][i] = A[i];
118                }
119                return d;
120        }
121
122        /**
123         *
124         * @param A
125         * @return
126         */
127        public static double[][] abs(double[][] A) {
128                int m = A.length;
129                int n = A[0].length;
130                double[][] C = new double[m][n];
131                for (int i = 0; i < m; i++) {
132                        for (int j = 0; j < n; j++) {
133                                C[i][j] = Math.abs(A[i][j]);
134                        }
135                }
136                return C;
137        }
138
139        /**
140         *
141         * @param A
142         * @return
143         */
144        public static double[] abs(double[] A) {
145                int m = A.length;
146
147                double[] C = new double[m];
148                for (int i = 0; i < m; i++) {
149                        C[i] = Math.abs(A[i]);
150                }
151                return C;
152        }
153
154        // return a random m-by-n matrix with values between 0 and 1
155        /**
156         *
157         * @param m
158         * @param n
159         * @return
160         */
161        public static double[][] random(int m, int n) {
162                double[][] C = new double[m][n];
163                for (int i = 0; i < m; i++) {
164                        for (int j = 0; j < n; j++) {
165                                C[i][j] = Math.random();
166                        }
167                }
168                return C;
169        }
170
171        // return n-by-n identity matrix I
172        /**
173         *
174         * @param n
175         * @return
176         */
177        public static double[][] identity(int n) {
178                double[][] I = new double[n][n];
179                for (int i = 0; i < n; i++) {
180                        I[i][i] = 1;
181                }
182                return I;
183        }
184
185        // return x^T y
186        /**
187         *
188         * @param x
189         * @param y
190         * @return
191         */
192        public static double dot(double[] x, double[] y) {
193                if (x.length != y.length) {
194                        throw new RuntimeException("Illegal vector dimensions.");
195                }
196                double sum = 0.0;
197                for (int i = 0; i < x.length; i++) {
198                        sum += x[i] * y[i];
199                }
200                return sum;
201        }
202
203        // return C = A^T
204        /**
205         *
206         * @param A
207         * @return
208         */
209        public static double[][] transpose(double[][] A) {
210                int m = A.length;
211                int n = A[0].length;
212                double[][] C = new double[n][m];
213                for (int i = 0; i < m; i++) {
214                        for (int j = 0; j < n; j++) {
215                                C[j][i] = A[i][j];
216                        }
217                }
218                return C;
219        }
220
221        // return C = A + B
222        /**
223         *
224         * @param A
225         * @param B
226         * @return
227         */
228        public static double[][] add(double[][] A, double[][] B) {
229                int m = A.length;
230                int n = A[0].length;
231                double[][] C = new double[m][n];
232                for (int i = 0; i < m; i++) {
233                        for (int j = 0; j < n; j++) {
234                                C[i][j] = A[i][j] + B[i][j];
235                        }
236                }
237                return C;
238        }
239
240        // return C = A - B
241        /**
242         *
243         * @param A
244         * @param B
245         * @return
246         */
247        public static double[][] subtract(double[][] A, double[][] B) {
248                int m = A.length;
249                int n = A[0].length;
250                double[][] C = new double[m][n];
251                for (int i = 0; i < m; i++) {
252                        for (int j = 0; j < n; j++) {
253                                C[i][j] = A[i][j] - B[i][j];
254                        }
255                }
256                return C;
257        }
258
259        // return C = A * B
260        /**
261         *
262         * @param A
263         * @param B
264         * @return
265         */
266        public static double[][] multiply(double[][] A, double[][] B) {
267                int mA = A.length;
268                int nA = A[0].length;
269                int mB = B.length;
270                int nB = B[0].length;
271                if (nA != mB) {
272                        throw new RuntimeException("Illegal matrix dimensions.");
273                }
274                double[][] C = new double[mA][nB];
275                for (int i = 0; i < mA; i++) {
276                        for (int j = 0; j < nB; j++) {
277                                for (int k = 0; k < nA; k++) {
278                                        C[i][j] += (A[i][k] * B[k][j]);
279                                }
280                        }
281                }
282                return C;
283        }
284
285        // matrix-vector multiplication (y = A * x)
286        /**
287         *
288         * @param A
289         * @param x
290         * @return
291         */
292        public static double[] multiply(double[][] A, double[] x) {
293                int m = A.length;
294                int n = A[0].length;
295                if (x.length != n) {
296                        throw new RuntimeException("Illegal matrix dimensions.");
297                }
298                double[] y = new double[m];
299                for (int i = 0; i < m; i++) {
300                        for (int j = 0; j < n; j++) {
301                                y[i] += (A[i][j] * x[j]);
302                        }
303                }
304                return y;
305        }
306
307        /**
308         *
309         * @param A
310         * @param x
311         * @return
312         */
313        public static double[][] scale(double[][] A, double[] x) {
314                int m = A.length;
315                int n = A[0].length;
316                if (x.length != A.length) {
317                        throw new RuntimeException("Illegal matrix dimensions.");
318                }
319                double[][] y = new double[m][n];
320                for (int i = 0; i < m; i++) {
321                        for (int j = 0; j < n; j++) {
322                                y[i][j] = A[i][j] * x[i];
323                        }
324                }
325                return y;
326        }
327
328        /**
329         *
330         * @param A
331         * @param x
332         * @return
333         */
334        public static double[][] scale(double[][] A, double x) {
335                int m = A.length;
336                int n = A[0].length;
337
338                double[][] y = new double[m][n];
339                for (int i = 0; i < m; i++) {
340                        for (int j = 0; j < n; j++) {
341                                y[i][j] = A[i][j] * x;
342                        }
343                }
344                return y;
345        }
346
347
348        // vector-matrix multiplication (y = x^T A)
349        /**
350         *
351         * @param x
352         * @param A
353         * @return
354         */
355        public static double[] multiply(double[] x, double[][] A) {
356                int m = A.length;
357                int n = A[0].length;
358                if (x.length != m) {
359                        throw new RuntimeException("Illegal matrix dimensions.");
360                }
361                double[] y = new double[n];
362                for (int j = 0; j < n; j++) {
363                        for (int i = 0; i < m; i++) {
364                                y[j] += (A[i][j] * x[i]);
365                        }
366                }
367                return y;
368        }
369
370        // test client
371        /**
372         *
373         * @param args
374         */
375        public static void main(String[] args) {
376                System.out.println("D");
377                System.out.println("--------------------");
378                double[][] d = {{1, 2, 3}, {4, 5, 6}, {9, 1, 3}};
379                StdArrayIO.print(d);
380                System.out.println();
381
382                System.out.println("I");
383                System.out.println("--------------------");
384                double[][] c = Matrix.identity(5);
385                StdArrayIO.print(c);
386                System.out.println();
387
388                System.out.println("A");
389                System.out.println("--------------------");
390                double[][] a = Matrix.random(5, 5);
391                StdArrayIO.print(a);
392                System.out.println();
393
394                System.out.println("A^T");
395                System.out.println("--------------------");
396                double[][] b = Matrix.transpose(a);
397                StdArrayIO.print(b);
398                System.out.println();
399
400                System.out.println("A + A^T");
401                System.out.println("--------------------");
402                double[][] e = Matrix.add(a, b);
403                StdArrayIO.print(e);
404                System.out.println();
405
406                System.out.println("A * A^T");
407                System.out.println("--------------------");
408                double[][] f = Matrix.multiply(a, b);
409                StdArrayIO.print(f);
410                System.out.println();
411        }
412}