MCR Fails When All-zero Spectral Channels Present

[AndrewHerzing]

@CCampJr : @jat255 and I have been testing our Hyperspy MCR implementation. We noticed that if one or more spectral channels in the dataset contains all zeros then MCR returns empty factors. Is this behavior to be expected? If so, is there an accepted way of working around the limitation?

[CCampJr]

@AndrewHerzing @jat255

In general, all 0's is a no-no (in all sorts of optimization methods) because there is no instability to cause the algorithm to move/evolve. Better is to seed with small random values.

I could foresee adding a check-for-constants/0's and give a warning to the log

[jat255]

Do you see any negative consequences of adding an infinitesimally-small value to channels that are uniformly-zero before proceeding with the MCR algorithm? I suppose it could preferentially bias the solution positive or negative?

[charlesll]

@jat255 why not just drop those channels? They are virtually information-less for the algorithm, so no need to include them I think.

[charlesll]

@CCampJr we could write an example on this subject!

[CCampJr]

@AndrewHerzing @jat255 I may have misinterpreted "spectral channel". Do you mean portions of spectra or do you mean, e.g., you seeded in 3 initial guess spectra and 2 of them were just 0's?

@charlesll I agree!

[jat255]

@CCampJr I mean portions of the spectra. i.e. if you have an x y spectrum image map of spectra, if one of the channels is zero-valued at every x and y position, then you'll get the error described.

@charlesll I don't think we can "drop" channels (at least in the HyperSpy sense). If for instance, you have a dead pixel on your detector or something, removing the channel entirely would break the axis calibration. Admittedly, having all zero values is probably a pretty rare occurrence, but it could feasibly happen in very sparse maps, or if your detector has some incorrect offset or something (which was the case that caused me to notice the problem).

I don't know MCR intimately, but my hunch is that the best solution will just to make the zero channels actually be infinitesimally small (but definite), so that the algorithm can proceed.

[CCampJr]

@jat255 So by channels you mean a particular energy bin on the CCD? So across the entire image, with (e.g.) 1024 energy/frequencies channels/bins, it is identically zero across the whole image?

Just checking because I'll test this real quick to see what's up. Making it a very small random value would work, but I'm a bit surprised at the problem now that I understand what's going on.

[CCampJr]

@jat255 So I did not have this problem. See the images below. I set the 50th and 100th spectral bin to be identically 0 at every image pixel. Note the "max spectrum..." image is the maximum spectral intensity across the entire HSI image (I'm using the Demo dataset in the Jupyter NB)

[jat255]

I just saw you responded, but I'll send this comment through anyway...

Yes, zero across the whole image. Here's an example of what I mean (plotted within HyperSpy). This is two signals made of three gaussians with added noise. In one of them, I set the 25th energy channel to zero. This is a "pixel-first" view of the data (i.e. a map of spectra):

And this is "energy-first" (i.e. a series of images):

When I run our "decomposition" routine using pyMCR, I get a valid result for the "no-zeros" case, but not for the one with one channel forced to zero:

[jat255]

I don't believe HyperSpy is doing anything unusual, but I will try to reproduce using only pyMCR and not HyperSpy to make sure I see the same behavior... Give me a few to work it out.

[CCampJr]

Here's a more intense example (I randomly killed half of the channels). Just for info I used:

# Just a copy of my hyperspectral image data
hsi2 = 1*hsi

# Kill half the spectral pixels across the image
locs_to_kill = np.random.randint(0,200,100)
hsi2[...,locs_to_kill] = 0

initial_spectra2 = 1*initial_spectra
initial_spectra2[...,locs_to_kill] = 0

mcrar = McrAR(max_iter=100, st_regr='OLS', c_regr='OLS', 
              c_constraints=[ConstraintNonneg(), ConstraintNorm()],
              st_constraints=[ConstraintNonneg()])

mcrar.fit(hsi2.reshape((-1, wn.size)), ST=initial_spectra2, verbose=True)
print('\nFinal MSE: {:.7e}'.format(mcrar.err[-1]))

[CCampJr]

@jat255

I don't believe HyperSpy is doing anything unusual, but I will try to reproduce using only pyMCR and not HyperSpy to make sure I see the same behavior... Give me a few to work it out.

Maybe it's a dumb question, but are your spectra/images integer type? Weird things happen if you put in an image with integer type and then work in float space.

[jat255]

Unfortunately, no they're floats. It's required by the SVD that we run prior to feeding those in as initial conditions to pyMCR

[AndrewHerzing]

We should clarify that in our implementation of MCR, we first perform an SVD decomposition using numpy.linalg.svd. The resulting factors are then rotated by varimax before being given as the input for PyMCR.

[CCampJr]

@AndrewHerzing @jat255

• How many SVD components are you keeping?
• Make sure the input spectra (into MCR) have positive values?
• If you all want, you can send me data (numpy binaries would be fine)
  • HSI
  • SVD-varimax'd initial spectra