# How to efficiently apply PCA followed by SVD to extract the components of PCA?

I am working on segmenting a textured image using Gabor Filters. This is based on this paper. There exist a Matlab implementation of this as well. Earlier I managed to do up to this. With the help of this now I come up with the following:

(* Design an array of Gabor Filters *)

{numRows, numCols} = ImageDimensions[img];

wavelengthMin = 4/Sqrt[2];
wavelengthMax = Sqrt[nrows^2 + ncols^2];
n = Floor[Log2[wavelengthMax/wavelengthMin]];
wavelength = Table[{2, 2} 2^i, {i, 0, n - 2}];
deltaTheta = Pi/4;
orientation = Table[i, {i, 0, Pi - deltaTheta, deltaTheta}];

gabormag =
Flatten[Table[
GaborFilter[img, 1, wavelength[[i]], orientation[[j]]]], {i, 1,
Length@wavelength}, {j, 1, Length@orientation}]];
gaborWavelength =
Flatten[ConstantArray[
Table[Norm[wavelength[[i]]], {i, 1, Length@wavelength}],
Length@orientation]];


(* Post-processing of Gabor Magnitude images into Gabor Features *)

K = 3;
gabormagfiltered = gabormag;
Table[
sigma = 0.5 gaborWavelength[[i]];
gabormagfiltered[[i]] =
ImageData[GaussianFilter[gabormag[[i]], K sigma]]
, {i, 1, Length@gabormag}];


(* Add a map of spatial location information *)

meshgrid[x_List, y_List] := {ConstantArray[x, Length[y]],
Transpose@ConstantArray[y, Length[x]]}

{X, Y} = meshgrid[Range[ncols], Range[nrows]];
featureset = Append[Append[gabormagfiltered, X], Y];


(* Reshape Data *)

npoints = nrows*ncols;
freq = Dimensions[featureset][[1]];

Needs["Developer"]
X = DeveloperToPackedArray[
ArrayReshape[featureset, {freq, npoints}]];


(* Normalize the features to be zero mean, unit variance *)

Xstandardized = Standardize[X];


(* The following is based on the answer provided by @Michael E2 here) ( Apply Principal Component Analysis *)

pca = PrincipalComponents[Xstandardized];
X1 = Transpose[Transpose@Xstandardized - Mean[Xstandardized]];


(* Extract coefficients *)

{u, \[Sigma], v} = SingularValueDecomposition[X1];


However, when I apply SingularValueDecomposition on X1 my machine (intel core i7, 16 GB RAM) runs out of memory, while in Matlab it works absolutely fine.

How can I efficiently apply PCA and SVD to extract the coefficients that I get directly using pca in Matlab?

For example, below are the original image, feature2DImage and clustered image (with 5 clusters)

• If you don't need all of them, try {u, \[Sigma], v} = SingularValueDecomposition[X1, Length@SingularValueList[X1]]; – Michael E2 May 8 '18 at 12:41
• @MichaelE2 After getting v following your approach, I applied feature2DImage = ImageAdjust@ImageResize[Image[X1.v], {225, 225}]. However, this does not give the desired image. Where am I doing wrong? – Majis May 8 '18 at 12:49
• @MichaelE2 Probably Dimensions[v] = {50625, 21} is a problem. – Majis May 8 '18 at 12:51
• I think what you want is PadRight[u.\[Sigma], {Automatic, Length@v}], since for the full SVD, u.\[Sigma] == X1.v. But it still doesn't look like the right image. If you have MATLAB, you might check the intermediate results. Somehow I doubt X1 is supposed to have a dimension 50625 if it's used to reconstruct the image as shown. -- Note Standardize zeros the mean, so Transpose[Transpose@Xstandardized - Mean[Xstandardized]] is unnecessary. – Michael E2 May 8 '18 at 13:34
• Try SingularValueDecomposition[X1, Min@Dimensions@X1] and u.\[Sigma] for coeff, then. -- It would help if there were an example image that everybody could do the same computations as each other. I don't think the web image of the puppy is exactly the same, or maybe you're using a different image. – Michael E2 May 8 '18 at 14:31