ICA  MATLAB  ASSIGNMENT

                               Michael Zibulevsky



*** Matlab functions are written in CAPITAL letters ***

1 - Qiasi-maximum likelihood source separation using the optimization toolbox function FMINUNC.m
--------------------------------------------------------------------------------------------

* Create N 1D sparse sources of length T using SPRANDN function.
  (N=2-10; T=500).

* Create N 1D uniformly distributed sources of length T with zero mean :
   
            S=rand(N,T) - 0.5;


* Create a mixing matrix using RAND, and create mixtures: X=AS.

* Define the function h(x) - "approximate log-pdf", both for
  supergaussian and subgaussian sources:

     supergaussian (for sparse signals)- smooth approximation of the absolute value 
                                            h(x)=SQRT(x.^2+lam^2), lam~0.1
     subgaussian (for uniform like sources) 
                                            h(x)=abs(x).^k, k~4

* Write a Matlab function that computes the minus-loglikelihood according to the
  paper "Blind Source Separation with Relative Newton Method"
  (this function  will be given as a parameter to FMINUNC).

* Use FMINUNC to find the minimum of minus-loglikelihood using only the function value.
  Find unmixing matrix W and estimated sources. The estimated sources rec=inv(W1)*X.

* Add the calculations of Gradient and Hessian to the function and test
  the separation results with only fun+grad and fun+grad+hess.
  separation results can be tested with:
  a) W*A - should have only 1 big component in every line.
  b) plot reconstructed sources over original sources.



* How to check that your Gradient and Hessian are correct:
 
    let f(x), g(x) and H(x)  be a function, its gradient and Hessian correspondingly.
    ( In our project x=W(:) )

    We can compute numerically the approximation of

        i-th element of gradient:  g_i ~ (f(x+ eps*e_i) -f(x))/eps
	i-th column of Hessian     H_i ~ (g(x+ eps*e_i) -g(x))/eps

    with eps=1e-7;

   Check that the numerical gradient and Hessian are very close to the analytical ones,
   within error less then 1e-5, for random values of x.
  

2 - Newton method with "frozen" Hessian 
---------------------------------------
* Download the newt_froz code from  http://ie.technion.ac.il/~mcib/
* Take images and create mixtures of them (same as above).
* Since the images are not usually sparse, we should make
  them sparse. The most simple way is to use GRADIENT and
  concatenate the derivatives in X and Y directions to one
  "mixture".
* Now we have sparse mixtures, we can use the function from
  assignment 1 as an input to NEWT_FROZ.
* Find the unmixing matrix using NEWT_FROZ and operate it on
  the ORIGINAL images (before gradient).
* Test the separation results (W*A or any other way).

3 - Relative Newton method 
---------------------------------------
* Download Relative Newton code from  http://ie.technion.ac.il/~mcib/
* Same as newt_froz. Compare the time to convergence..

BONUS 1 - Optimal basis for sounds
---------------------------------------
* Resample a wav file in matlab using RESAMPLE to 8000Hz.
* Do the next for K=[64,128] R=[1^4,10^5] (all combinations):
* Take a part of the wav in length K*R samples and reshape to
  the matrix X of size RxK, forming R samples of length K.
* Find the optimal (most sparse) basis using RELNEWTON.
* Examine the histogram of X vs. the histogram of the coefficients C.
* Examine the imagesc(abs(X)) vs.  imagesc(abs(C)).
* See that the basis PHI contains harmonic functions (sinus in various
  frequencies).

BONUS 2 - Comparison between relnewton and fastICA
---------------------------------------
* For the sound separation see if fastICA converges quicker.
  Compare the error for both algorithms.
* Test if using fastICA as an initial solution improves the convergence
  time.

GOOD LUCK!!