CN105138826A - Raman signal reconstruction method under strong noise background - Google Patents

Abstract

The invention discloses a raman signal reconstruction method under a strong noise background. From the perspective of signal reconstruction, a useful raman signal is reconstructed from an original measurement signal under the strong noise background by gaussian analysis dictionary-based signal sparse representation; noise removal and baseline correction can be simultaneously achieved by the reconstructed signal; meanwhile, through calculating the inner product of the dictionary and the residual error, quickening of the signal processing rate and reduction of memory usage can be achieved; effective reconstruction of the raman signal under the strong noise background can be achieved by the algorithm; a raman peak is accurate to displace and position; the peak area attenuation is small; the signal to nose ratio is obviously improved; the algorithm processing efficiency is high; and the raman signal reconstruction method can be applied to signal preprocessing of a high temporal-spatial resolution raman measurement system.

Description

Raman signal reconstructing method under a kind of strong noise background

Technical field

The present invention relates to signal transacting field, the Raman signal reconstructing method particularly under a kind of strong noise background.

Background technology

Raman spectroscopy is a kind of strong molecular structure research tool, but Raman scattering effect is very faint, and its scattered light intensity is about 10 of incident intensity ^-8-10 ^-6.In order to obtain high-spatial and temporal resolution, under usually adopting very little integral time (as 0.1s) situation, now, under faint effective Raman signal can be submerged in strong noise background.

Conventional signal antinoise method has moving window average smoothing method, EMD (EmpiricalModeDecomposition empirical mode decomposition), Wavelet Transform etc.

Moving window average smoothing method, based on the difference between signal and noise statistics, its basic assumption is noise is zero mean noise, and reaching the object improving signal to noise ratio (S/N ratio) by averaging to original signal, is stress release treatment most popular method.Be the smooth window of odd number 2 ω+1 by selecting a width, with centre wavelength point k for reference point from left to right moving window, the mean value of measurements all in forms overlay area is replaced the measured value that centre wavelength point is corresponding, until complete a little level and smooth, the method Window width ω affects sharpening result, and width is too little, smooth effect is not good, width is too large, then smooth out characteristic peak information, cause spectrum distortion; And there is boundary problem.For the smooth window that width is 2 ω+1, two ends, spectrum left and right respectively have ω point not to be processed;

EMD by signal adaptive resolve into limited from high to low, ascending IMF (instrinsicmodefunction intrinsic mode function) component of time scale and single trend term, but lack strict mathematics to be difficult to determine according to, stopping criterion, and there is end effect, modal overlap effect, algorithm stability is bad.

Wavelet Transform is via the flexible of wavelet basis and translation, reach the localization of signal time frequency analysis, can the simultaneously temporal signatures of stick signal and frequency domain character, in suitable yardstick, the effective constituent of non-stationary signal can expose and noise different characteristics, uses signal and the various transmission characteristics of noise in multiple dimensioned scope, useful signal can be gone out in noise background extracting, but choosing denoising result interference of wavelet basis is comparatively large, and simultaneously when low signal-to-noise ratio, effect is poor.

During traditional denoising method process low signal-to-noise ratio Raman signals such as moving window average smoothing method, EMD, Wavelet Transform, effect is unsatisfactory, can not realize effective reduction of Raman signal under strong background noise.

Summary of the invention

Fundamental purpose of the present invention is to overcome above-mentioned defect of the prior art, proposes the Raman signal reconstructing method under a kind of strong noise background, and its algorithm process efficiency is high, signal noise ratio improve realizes noise remove and baseline correction obviously, simultaneously.

The present invention adopts following technical scheme:

A Raman signal reconstructing method under strong noise background, is characterized in that, pre-defined initial auxiliary inner product is α ⁰=D ^ty, y are original Raman Measurement signal, and D is Gauss's dictionary; Companion matrix G=D ^td; Initial square error ε ⁰=y ^ty; Target square error is ε; Residual error r=(y-D γ), γ is the sparse expression of y; Comprise the steps:

1) initialization: make iterative process selected atom sequence I:=(), to companion matrix G _i,Icholesky decompose L initialize L:=[1], auxiliary inner product α :=α ⁰, iterative process intermediate variable δ ⁰:=0, primary iteration sequence number n:=1;

2) calculate the auxiliary inner product α of Gauss's dictionary atom and residual signals, select the atom index mated most with current residue signal judge whether current iteration sequence number is greater than 1, if so, then enters step 3); Otherwise, jump to step 4);

3) L is upgraded by introducing intermediate variable;

4) the atom index mated most will selected add set

5) according to formula for D _i ^ty is the submatrix D of Gauss's dictionary D of the correspondence of the atomic series composition of current renewal _iinner product with tested original Raman signal y, solves c, c be upgrade after current Gauss's dictionary to the sparse expression of original Raman Measurement signal, i.e. γ _i;

6) the iteration intermediate variable β=G after n-th iteration is calculated according to the atomic series after renewal _iγ _i, wherein G _ifor the atomic series companion matrix G=D after correspondence renewal ^tthe submatrix of D, and make current auxiliary inner product α equal initially to assist inner product α ⁰deduct the iteration intermediate variable β obtained;

7) another intermediate variable after n-th iteration is tried to achieve according to the atomic series after renewal β _ithe subvector of the iteration intermediate variable β corresponding to the atomic series of renewal, upgrades n-th iteration square error ε ⁿ=ε ^n-1-δ ⁿ+ δ ^n-1;

8) judge whether to meet end condition, if so, then export original Raman Measurement signal y sparse expression γ, effective Raman signal x=D γ of reconstruct, otherwise n=n+1, return step 2).

Preferably, described step 3) be specially, be that k atom adds last row of Gauss's dictionary D to by the sequence number of trying to achieve, upgrade companion matrix basis again ask for intermediate variable ω, according to the ω asked for and formula $L:=\left[\begin{array}{cc}L& 0\\ {\mathrm{\ω}}^T& \sqrt{1-{\mathrm{\ω}}^T\mathrm{\ω}}\end{array}\right],$ Upgrade L.

Preferably, in step 8) in, described end condition is the square error ε of n-th iteration ⁿ< ε.

From the above-mentioned description of this invention, compared with prior art, the present invention has following beneficial effect:

The angle that the present invention is based on signal reconstruction extracts effective Raman signal from raw measured signal, the curve of spectrum of reconstruct is smooth, can realize noise remove and baseline correction simultaneously, and Raman peaks displacement location is accurate, peak area decay is little, signal noise ratio improve is obvious, and algorithm process efficiency is high.

Accompanying drawing explanation

Fig. 1 is process flow diagram of the present invention;

Fig. 2 is desirable simulate signal;

Fig. 3 adds the signal after making an uproar to the desirable simulate signal of Fig. 2;

Fig. 4 is the contrast of emulation noisy signal after Batch-OMP process and desirable simulate signal;

Fig. 5 is the effect contrast figure of the emulation adopting moving window average, wavelet transformation and the original noisy signal of Batch-OMP process respectively;

Fig. 6 is the effect contrast figure of the emulation not adopting moving window average, wavelet transformation and Batch-OMP process original signal.

Embodiment

Below by way of embodiment, the invention will be further described.

Raman signal reconstructing method under a kind of strong noise background of the present invention, its principle adopted picks out the maximum candidate's atom of inner product by often walking, again selected atom is carried out Schimidt orthogonalization process, signal being projected in these orthogonal atomic forms spatially, obtain the component of signal on candidate's atom and remnants, thereafter adopt identical method to decompose residual component, entered limited number of time and decompose, and obtained original signal former linearly molecular by limited.The step of Schimidt orthogonalization process ensure that the optimality of single iteration, but the iterative process of least square method is comparatively complicated, and therefore OMP algorithm exists the problem that solving speed is comparatively slow and memory space is larger.The present invention decomposes by Cholesky and calculating inner product assists inner product D ^tr replaces directly calculating residual error r, and improve algorithm process efficiency, the concrete implementing procedure block diagram of Batch-OMP of the present invention as shown in Figure 1.Sparse expression dictionary can be divided into the analysis dictionary based on specific function and the study dictionary based on sample learning.Under useful Raman signal is submerged in strong background noise, according to sample learning dictionary, noise can be expressed as useful signal by mistake, therefore our Analysis about Selection dictionary.Because Raman peaks can be approximately Gaussian peak, and each parameter of Gaussian function (peak intensity, peak position and standard deviation) explicit physical meaning, therefore we select Gaussian function to carry out sparse expression as analysis dictionary to the Raman signal under strong noise background.

At the end of defining n-th iteration, the sparse expression of residual signals and signal y is respectively r ⁿ, γ ⁿ, Gauss's dictionary is D, because each iteration residual error is all orthogonal with current selected atom, so for all iteration sequence number n, meets (r ⁿ) ^td γ ⁿ=0

r ⁿ＝y-Dγ ⁿ＝y-Dγ ^n-1+Dγ ^n-1-Dγ ⁿ＝r ^n-1+D(γ ^n-1-γ ⁿ)

Definition square error and introduce companion matrix G=D ^td, can be obtained fom the above equation

$\begin{array}{c}\left|\right|r^n|{|}_2^2={\left(r^n\right)}^T{r}^n\\ ={\left(r^n\right)}^T(r^{n-1}+D({\mathrm{\&gamma;}}^{n-1}-{\mathrm{\&gamma;}}^n\left)\right)\\ {\mathrm{=(r}}^{n-1}\mathrm{+D}\left({\mathrm{\&gamma;}}^{n-1}-{\mathrm{\&gamma;}}^n\right){)}^T{r}^{n-1}+{\left(r^n\right)}^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-{\left(r^{n-1}\right)}^T{\mathrm{D\&gamma;}}^n+{\left(r^n\right)}^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-{(y-{\mathrm{D\&gamma;}}^{n-1})}^T{\mathrm{D\&gamma;}}^n+{(y-{\mathrm{D\&gamma;}}^{n-1})}^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-y^T{\mathrm{D\&gamma;}}^n+y^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-{(r^n+{\mathrm{D\&gamma;}}^n)}^T{\mathrm{D\&gamma;}}^n+{(r^{n-1}+{\mathrm{D\&gamma;}}^{n-1})}^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-{\left({\mathrm{\&gamma;}}^n\right)}^T{D}^T{\mathrm{D\&gamma;}}^n+{\left({\mathrm{\&gamma;}}^{n-1}\right)}^T{D}^T{\mathrm{D\&gamma;}}^{n-1}\\ =\left|\right|r^{n-1}|{|}_2^2-{\left({\mathrm{\&gamma;}}^n\right)}^T{\mathrm{G\&gamma;}}^n+{\left({\mathrm{\&gamma;}}^{n-1}\right)}^T{\mathrm{G\&gamma;}}^{n-1}\end{array}$

Introduce intermediate variable δ ⁿ=(γ ⁿ) ^tg γ ⁿ, the square error step of updating therefore after n-th iteration can be expressed as ε ⁿ=ε ^n-1-δ ⁿ+ δ ^n-1, definition intermediate variable β=G γ; Due to the sparse expression γ=D of signal y ⁺y=(D ^td) ^-1d ^toriginal Raman Measurement signal sparse expression γ after y, introducing Cholesky decomposition solves n-th iteration _i, and then improve algorithm counting yield.

Of the present invention as follows about Batch-OMP algorithm idiographic flow:

Pre-defined: initial alpha ⁰=D ^ty, y are original Raman Measurement signal, and D is Gauss's dictionary, companion matrix G=D ^td, initial square error ε ⁰=y ^ty, target square error ε.

1) initialization: iterative process selected atom sequence I:=(), companion matrix G _i,Icholesky decompose initialize L:=[1], auxiliary inner product α :=α ⁰, iterative process intermediate variable δ ⁰:=0, primary iteration sequence number n:=1;

2) the auxiliary inner product α=D of each atom of Dictionary of Computing and residual signals ^tr, the absolute value of the atom of trying to achieve and residual signals inner product of vectors is larger, represent that selected atom more mates with residual error structure in current iteration computing, so select the atom maximum with current iteration residual error inner product to carry out projecting and just can realize farthest approaching current residue signal with minimum nonzero coefficient, therefore tentatively select the atom index mated most with residual signals judge whether current iteration sequence number is greater than 1, if current iteration sequence number is greater than 1, then jump to step 3); Otherwise, jump to step 4).

3) by step 2) sequence number of trying to achieve is that k atom adds last row of dictionary D to, upgrades companion matrix basis again ask for intermediate variable ω, then according to formula $L:=\left[\begin{array}{cc}L& 0\\ {\mathrm{\&omega;}}^T& \sqrt{1-{\mathrm{\&omega;}}^T\mathrm{\&omega;}}\end{array}\right]$ (L is G to upgrade L _i,Icholesky decompose, i.e. G _i,I=LL ^t, G _i,Ifor the submatrix of companion matrix G corresponding to current atomic series, this step is to reduce step 5) calculated amount and the space hold intermediate variable that calculates and introduce).

4) the atom index mated most with residual signals (i.e. original Raman Measurement signal y) before iteration first will be selected $\hat{k}:=\arg \underset{k}{max}\left\{\right|{\mathrm{\&alpha;}}_k\left|\right\}$ Add set $I:=(I,\hat{k}).$

5) after upgrading iteration atomic series, according to formula for D _i ^ty (the submatrix D of the dictionary D of the correspondence of the atomic series composition of current renewal _iwith the inner product of tested original Raman signal y), the sparse expression γ of the original Raman Measurement signal of current dictionary sequence pair is upgraded by solving c _i, c be upgrade after current Gauss's dictionary to the sparse expression of original Raman Measurement signal, i.e. γ _i.

6) adopt the atomic series after upgrading, calculate the iteration intermediate variable β=G after n-th iteration _iγ _i, wherein G _ifor the atomic series companion matrix G=D after correspondence renewal ^tthe submatrix of D, with initial auxiliary inner product α ⁰deduct iteration intermediate variable β and upgrade current auxiliary inner product α, i.e. α :=α ⁰-β.

7) n-th iteration intermediate variable is tried to achieve by the atom dictionary sequence after renewal β _ithe subvector of the iteration intermediate variable β corresponding to the atomic series of renewal; Upgrade n-th iteration square error ε ⁿ=ε ^n-1-δ ⁿ+ δ ^n-1.

8) if meet end condition: the square error ε that n-th iteration terminates ⁿ< ε (target square error), then exit circulation, exports original Raman Measurement signal y sparse expression γ, effective Raman signal x=D γ of reconstruct; Otherwise n=n+1, returns step 2), enter next iteration.

To in confirmatory experiment of the present invention, adopt normal signal denoising method moving window average smoothing method, Noise Elimination from Wavelet Transform method and employing of the present invention based on Batch-OMP denoise algorithm equal signal to noise ratio (S/N ratio) emulation Raman signal, test Raman signal and compare, result is as shown in Fig. 4, Fig. 5, Fig. 6, table 1, table 2, table 3.Fig. 2 is desirable simulate signal, and Fig. 3 adds the signal after making an uproar to the desirable simulate signal of Fig. 2, can find out that ideal signal is submerged in noise completely.Table 1 adds to emulation the contrast that Raman signal of making an uproar carries out the Raman shift after denoising and desirable simulate signal for Batch-OMP, can know and find out, the Raman signal peak shift accurate positioning after Batch-OMP denoising.Table 2, table 3 are moving window average, wavelet transformation and Batch-OMP process simulate signal signal to noise ratio snr and mean square deviation (MSE) improve comparing result, can find out, the inventive method (Batch-OMP) is compared with moving window average smoothing method, wavelet transformation, and signal to noise ratio (S/N ratio) (SNR) and mean square deviation (MSE) are improved obviously.Fig. 4 is the contrast of emulation noisy signal after Batch-OMP process and desirable simulate signal, can obviously find out, this method can reduce effective Raman signal, and peak area decay is little, can effectively reconstruct useful Raman signal.Fig. 5, Fig. 6 are respectively the effect contrast figure of moving window average, wavelet transformation and Batch-OMP process emulation noisy signal, experimental signal, we can find out, the inventive method relatively and additive method, denoising effect is better, and the signal adopting this method (Batch-OMP) to reconstruct can realize noise remove and baseline correction simultaneously.And the inventive method structure in specific implementation is simple, relative to conventional algorithm, does not increase too many complexity.

Table 1

	200cm -1	500cm -1	800cm -1
				Ideal signal	200	500	800
Batch-OMP	200	499	801

Table 2

	SNR (signal to noise ratio (S/N ratio) dB)	MSE (mean square deviation)
			Originally add simulate signal of making an uproar	-24.308	4.611e+03
Moving window average	-24.176	4.542e+03
			Wavelet transformation	-24.114	4.509e+03

Batch-OMP

5.779

0.144e+03

Table 3

	SNR (signal to noise ratio (S/N ratio) dB)	Mean square deviation MSE
			Initial experiments signal	-5.5019	4.2842
Moving window average	-1.5126	2.7036
			Wavelet transformation	2.3298	1.7371
Batch-OMP	4.7163	1.3198

The present invention also can be used for the Raman image of reduction low signal-to-noise ratio.Raman signal denoising shown in embodiment only reference as an example.For raising Raman image resolution, by peak position parameter and the restriction degree of rarefication of adjustment Gauss dictionary atom, also can realize reducing to the Raman image of low signal-to-noise ratio.

Above are only the specific embodiment of the present invention, but design concept of the present invention is not limited thereto, all changes utilizing this design the present invention to be carried out to unsubstantiality, all should belong to the behavior of invading scope.

Claims (3)

1. the Raman signal reconstructing method under strong noise background, is characterized in that, pre-defined initial auxiliary inner product is α ⁰=D ^ty, y are original Raman Measurement signal, and D is Gauss's dictionary; Companion matrix G=D ^td; Initial square error ε ⁰=y ^ty; Target square error is ε; Residual error r=(y-D γ), γ is the sparse expression of y; Comprise the steps:

1) initialization: make iterative process selected atom sequence I:=(), to companion matrix G _i,Icholesky decompose L initialize L:=[1], auxiliary inner product α :=α ⁰, iterative process intermediate variable δ ⁰:=0, primary iteration sequence number n:=1;

2) calculate the auxiliary inner product α of Gauss's dictionary atom and residual signals, select the atom index mated most with current residue signal judge whether current iteration sequence number is greater than 1, if so, then enters step 3); Otherwise, jump to step 4);

3) L is upgraded by introducing intermediate variable;

4) the atom index mated most will selected add set

5) according to formula for D _i ^ty is the submatrix D of Gauss's dictionary D of the correspondence of the atomic series composition of current renewal _iinner product with tested original Raman signal y, solves c, c be upgrade after current Gauss's dictionary to the sparse expression of original Raman Measurement signal, i.e. γ _i;

6) the iteration intermediate variable β=G after n-th iteration is calculated according to the atomic series after renewal _iγ _i, wherein G _ifor the atomic series companion matrix G=D after correspondence renewal ^tthe submatrix of D, and make current auxiliary inner product α equal initially to assist inner product α ⁰deduct the iteration intermediate variable β obtained;

7) another intermediate variable after n-th iteration is tried to achieve according to the atomic series after renewal β _ithe subvector of the iteration intermediate variable β corresponding to the atomic series of renewal, upgrades n-th iteration square error ε ⁿ=ε ^n-1-δ ⁿ+ δ ^n-1;

8) judge whether to meet end condition, if so, then export original Raman Measurement signal y sparse expression γ, effective Raman signal x=D γ of reconstruct, otherwise n=n+1, return step 2).

2. the Raman signal reconstructing method under a kind of strong noise background as claimed in claim 1, is characterized in that: described step 3) be specially, be that k atom adds last row of Gauss's dictionary D to by the sequence number of trying to achieve, upgrade companion matrix basis again ask for intermediate variable ω, according to the ω asked for and formula $L:=\left[\begin{array}{cc}L& 0\\ {\mathrm{\&omega;}}^T& \sqrt{1-{\mathrm{\&omega;}}^T\mathrm{\&omega;}}\end{array}\right],$ Upgrade L.

3. the Raman signal reconstructing method under a kind of strong noise background as claimed in claim 1, is characterized in that: in step 8) in, described end condition is the square error ε of n-th iteration ⁿ< ε.