Java Singular.fullSVD方法代码示例

import org.jblas.Singular; //导入方法依赖的package包/类
public Matrix[] calc(Matrix source) {
	final DoubleMatrix matrix;
	if (source instanceof JBlasDenseDoubleMatrix2D) {
		matrix = ((JBlasDenseDoubleMatrix2D) source).getWrappedObject();
	} else if (source instanceof HasColumnMajorDoubleArray1D) {
		matrix = new JBlasDenseDoubleMatrix2D(MathUtil.longToInt(source.getRowCount()), MathUtil.longToInt(source
				.getColumnCount()), ((HasColumnMajorDoubleArray1D) source).getColumnMajorDoubleArray1D())
				.getWrappedObject();
	} else {
		matrix = new JBlasDenseDoubleMatrix2D(source).getWrappedObject();
	}
	final DoubleMatrix[] svd = Singular.fullSVD(matrix);
	final Matrix u = new JBlasDenseDoubleMatrix2D(svd[0]);
	DoubleMatrix sVector = svd[1];
	final Matrix s = new JBlasDenseDoubleMatrix2D(MathUtil.longToInt(source.getRowCount()),
			MathUtil.longToInt(source.getColumnCount()));
	for (int i = 0; i < sVector.rows; i++) {
		s.setAsDouble(sVector.get(i), i, i);
	}
	final Matrix v = new JBlasDenseDoubleMatrix2D(svd[2]);
	return new Matrix[] { u, s, v };
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
/**
 * 
 * 
 * @param A
 * @param dim
 * @return
 */
public DoubleMatrix pca(DoubleMatrix A, int dim) {
	
	logger.info("Type 'pca' started (target dimensions = "+dim+").");
	
	A = centerData(A);
	save("A", A);
	
	logger.info("Computing SVD...");
	DoubleMatrix[] usv = Singular.fullSVD(A);
	DoubleMatrix U = usv[0];
	DoubleMatrix S = usv[1];
	
	save("U", U);
	save("S", S);
	
	// 
	logger.info("Reducing U to Uk...");
	DoubleMatrix Uk = new DoubleMatrix(U.rows, dim);
	for(int i=0; i<dim; i++)
		Uk.putColumn(i, U.getColumn(i));
	save("Uk", Uk);
	
	// build S matrix
	logger.info("Reducing S to Sk...");
	DoubleMatrix Sk = new DoubleMatrix(dim, dim);
	for (int i = 0; i < dim; i++) {
		Sk.put(i, i, S.get(i));
	}
	save("Sk", Sk);
	
	// calculate principal component matrix...
	logger.info("Computing principal components...");
	DoubleMatrix B = Uk.mmul(Sk);
	save("B", B);
	
	return B;
	
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
/**
 * @param A
 * @param k
 * @return
 */
public DoubleMatrix compress(DoubleMatrix A, int k) {
	
	logger.info("Type 'compress' started...");

	A = centerData(A);
	
	DoubleMatrix[] usv = Singular.fullSVD(A);
	// n x n
	DoubleMatrix U = usv[0];
	// n x p
	DoubleMatrix S = usv[1];
	// p x p (straight)
	DoubleMatrix V = usv[2];
			
	// k x k
	DoubleMatrix Sk = new DoubleMatrix(U.columns, V.columns);
	for (int i = 0; i < k; i++) {
		Sk.put(i, i, S.get(i));
	}
	
	// 
	DoubleMatrix Aapprox = U.mmul(Sk).mmul(V.transpose());
	save("Aapprox", Aapprox);
	
	return Aapprox;
	
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
public DoubleMatrix reconstruct2(DoubleMatrix A, int k) {
		
		logger.info("Type 'reconstruct2' started...");

		A = centerData(A);
		save("A", A);
		
		DoubleMatrix[] usv = Singular.fullSVD(A);
		// n x n
		DoubleMatrix U = usv[0];
		// n x p
		DoubleMatrix S = usv[1];
		// p x p (straight)
		DoubleMatrix V = usv[2];
		
		save("U", U);
		save("S", S);
		save("V", V);
		
		for(int i = k; i < U.columns; i++)
			U.putColumn(i, DoubleMatrix.zeros(U.rows, 1));
		save("Uk", U);
		
		DoubleMatrix Sm = new DoubleMatrix(A.rows, A.columns);
		for (int i = 0; i < k; i++) {
			Sm.put(i, i, S.get(i));
		}
		save("Sk", Sm);
		
//		for (int i = k; i < V.rows; i++) {
//			V.putColumn(i, DoubleMatrix.zeros(1, V.columns));
//		}
//		visual("Vk", V);
 		
		// 
		DoubleMatrix Aapprox = U.mmul(Sm).mmul(V.transpose());
		save("Aapprox", Aapprox);
				
		return Aapprox;
	}
 

import org.jblas.Singular; //导入方法依赖的package包/类
/**
 * 初始化PCA
 * 
 * @param x 特征矩阵
 */
public PCA(FloatMatrix x) {
	super();
	this.x = x;

	this.m = x.rows;
	FloatMatrix Sigma = x.transpose().mmul(x).div(m);
	FloatMatrix[] result = Singular.fullSVD(Sigma);

	this.u = result[0].neg();// 与matlab保持一致
	this.s = result[1];
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
public static DoubleMatrix pinv(final DoubleMatrix A) {
  final DoubleMatrix[] usv = Singular.fullSVD(A);
  final DoubleMatrix s_inv = new DoubleMatrix(usv[1].length);
  for (int i = 0; i < usv[1].length; i++) {
    s_inv.put(i, Math.abs(usv[1].get(i)) < TOLERANCE ? 0.0 : 1.0 / usv[1].get(i));
  }
  return usv[2].mmul(DoubleMatrix.diag(s_inv)).mmul(usv[0].transpose());
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
public static DoubleMatrix[] svd(DoubleMatrix src)
{
    DoubleMatrix[] usv = Singular.fullSVD(src);
       DoubleMatrix S = usv[1];
       usv[1] = convertDiagMatrix(S, src.rows, src.columns);
       return usv;
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
@Test
public void testSVD()
{
	double[][] testData = { 
			{ 36, 49, 47, 11 }, 
			{ 2, 68, 27, 42 }, 
			{ 42, 25, 38, 3 } 
		};
	RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
	SingularValueDecomposition svd = new SingularValueDecomposition(matrix);
       System.out.println(svd.getU());
       System.out.println(svd.getS());
       System.out.println(svd.getV());
	
	DoubleMatrix[] usv = Singular.fullSVD(new DoubleMatrix(testData));
       System.out.println(usv[0]);
       System.out.println(usv[1]);
       System.out.println(usv[2]);

       DoubleMatrix U = usv[0];
       DoubleMatrix S = usv[1];
       DoubleMatrix V = usv[2];
       DoubleMatrix mt = new DoubleMatrix(3, 4);
       for (int i = 0; i < S.length; i++)
       {
           mt.put(i, i, S.get(i));
       }
       System.out.println(mt.toString().replace(";", "\n"));
       DoubleMatrix src = U.mmul(mt).mmul(V.transpose());
       System.out.println(src.toString().replace(";", "\n"));
       mt = Solve.pinv(mt);
       System.out.println(mt.toString().replace(";", "\n"));
}
 

import org.jblas.Singular; //导入方法依赖的package包/类
/**
 * @param A
 * @param k
 * @return
 */
public DoubleMatrix reconstruct(DoubleMatrix A, int k) {
	
	logger.info("Type 'reconstruct' started...");

	A = centerData(A);
	save("A", A);
	
	DoubleMatrix[] usv = Singular.fullSVD(A);
	// n x n
	DoubleMatrix U = usv[0];
	// n x p
	DoubleMatrix S = usv[1];
	// p x p (straight)
	DoubleMatrix V = usv[2];
	
	save("U", U);
	save("S", S);
	save("V", V);
	
	// n x k
	DoubleMatrix Uk = new DoubleMatrix(U.rows, k);
	for(int i=0; i<k; i++)
		Uk.putColumn(i, U.getColumn(i));
	save("Uk", Uk);
	
	// k x k
	DoubleMatrix Sk = new DoubleMatrix(k, k);
	for (int i = 0; i < k; i++) {
		Sk.put(i, i, S.get(i));
	}
	save("Sk", Sk);
	
	// p x k (straight)
	DoubleMatrix Vk = new DoubleMatrix(A.columns, k);
	for (int i = 0; i < k; i++) {
		Vk.putColumn(i, V.getColumn(i));
	}
	save("Vk", Vk);
		
	// 
	DoubleMatrix Aapprox = Uk.mmul(Sk).mmul(Vk.transpose());
	save("Aapprox", Aapprox);
	
	DoubleMatrix Areduced = new DoubleMatrix(A.rows, k);
	for(int i=0; i<k; i++)
		Areduced.putColumn(i, Aapprox.getColumn(i));
	save("C"+k, Areduced);
	
	return Areduced;
	
}