CN117574064A - Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method - Google Patents
Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method Download PDFInfo
- Publication number
- CN117574064A CN117574064A CN202311353099.9A CN202311353099A CN117574064A CN 117574064 A CN117574064 A CN 117574064A CN 202311353099 A CN202311353099 A CN 202311353099A CN 117574064 A CN117574064 A CN 117574064A
- Authority
- CN
- China
- Prior art keywords
- peak
- value
- raman spectrum
- axis
- follows
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000001237 Raman spectrum Methods 0.000 title claims abstract description 81
- 238000004364 calculation method Methods 0.000 title claims abstract description 21
- 238000000034 method Methods 0.000 claims abstract description 55
- 238000001228 spectrum Methods 0.000 claims abstract description 9
- 230000014509 gene expression Effects 0.000 claims description 30
- 230000003595 spectral effect Effects 0.000 claims description 26
- 238000009499 grossing Methods 0.000 claims description 17
- 238000001069 Raman spectroscopy Methods 0.000 claims description 4
- 239000011159 matrix material Substances 0.000 claims description 4
- 238000004451 qualitative analysis Methods 0.000 abstract description 2
- 238000004445 quantitative analysis Methods 0.000 abstract description 2
- 238000010586 diagram Methods 0.000 description 5
- 238000000926 separation method Methods 0.000 description 2
- 229920002430 Fibre-reinforced plastic Polymers 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 239000011151 fibre-reinforced plastic Substances 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
- G01N21/62—Systems in which the material investigated is excited whereby it emits light or causes a change in wavelength of the incident light
- G01N21/63—Systems in which the material investigated is excited whereby it emits light or causes a change in wavelength of the incident light optically excited
- G01N21/65—Raman scattering
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/10—Pre-processing; Data cleansing
- G06F18/15—Statistical pre-processing, e.g. techniques for normalisation or restoring missing data
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/08—Feature extraction
- G06F2218/10—Feature extraction by analysing the shape of a waveform, e.g. extracting parameters relating to peaks
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Physics (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Analysis (AREA)
- Computational Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Computation (AREA)
- Databases & Information Systems (AREA)
- Artificial Intelligence (AREA)
- Software Systems (AREA)
- Algebra (AREA)
- Probability & Statistics with Applications (AREA)
- Health & Medical Sciences (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Operations Research (AREA)
- Biochemistry (AREA)
- Analytical Chemistry (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Chemical & Material Sciences (AREA)
- General Health & Medical Sciences (AREA)
- Computing Systems (AREA)
- Investigating, Analyzing Materials By Fluorescence Or Luminescence (AREA)
Abstract
The invention discloses a Raman spectrum peak-splitting intensity calculation method based on Gaussian fitting, which comprises the following steps: step 1, acquiring an original spectrum image; step 2, acquiring a Raman spectrum curve; step 3, acquiring a Raman spectrum smooth curve; step 4, automatically obtaining peak dividing positions; step 5, acquiring an x-axis physical value of a Raman spectrum curve; step 6, initializing peak-splitting Gaussian parameters; step 7, fine tuning the split peak Gaussian parameters; and 8, acquiring the peak splitting area of the Raman spectrum. The method is suitable for the situation that personnel are difficult to interfere when calculating the peak splitting intensity of the Raman spectrum, and the accuracy of the automatically calculated peak splitting area of the Raman spectrum is higher; the method is also suitable for the condition of calculation of the peak intensity of the Raman spectrum, and has great application value for pure qualitative analysis, high-degree quantitative analysis and determination of molecular structures.
Description
Technical Field
The invention belongs to the technical field of Raman spectrum peak splitting intensity, and relates to a Raman spectrum peak splitting intensity calculation method based on Gaussian fitting.
Background
Although the raman spectroscopy technology has been widely used in various fields such as chemistry, physics, biology and medicine, the area calculation of raman spectroscopy peak separation still has some problems. Because of certain pollution of a target medium or an object, certain interference exists in the acquired Raman spectrum, so that certain difficulty exists in calculation of the peak splitting area of the Raman spectrum, manual intervention is generally needed, the accuracy and efficiency of a Raman spectrum peak splitting intensity calculation method are affected, and the problem of the application field of a computer is solved.
Disclosure of Invention
The invention aims to provide a Raman spectrum peak-splitting intensity calculation method based on Gaussian fitting, which solves the problems that in the area calculation process of Raman spectrum peak-splitting in the prior art, manual intervention is needed, the accuracy of a calculation result is not enough, and the calculation process is time-consuming.
The technical scheme adopted by the invention is that the Raman spectrum peak splitting intensity calculation method based on Gaussian fitting is implemented according to the following steps:
step 1, acquiring an original spectrum image;
step 2, acquiring a Raman spectrum curve;
step 3, acquiring a Raman spectrum smooth curve;
step 4, automatically obtaining peak dividing positions;
step 5, acquiring an x-axis physical value of a Raman spectrum curve;
step 6, initializing peak-splitting Gaussian parameters;
step 7, fine tuning the split peak Gaussian parameters;
and 8, acquiring the peak splitting area of the Raman spectrum.
The method has the beneficial effects that the method is suitable for the situation that personnel are difficult to intervene when calculating the peak splitting intensity of the Raman spectrum, and the accuracy of the automatically calculated peak splitting area of the Raman spectrum is higher; the method is also suitable for the situation of calculation of the peak intensity of the Raman spectrum, and has great application value for pure qualitative analysis, high-degree quantitative analysis and determination of molecular structures.
Drawings
FIG. 1 is a schematic overall flow diagram of the method of the present invention;
FIG. 2 is a schematic representation of a spectral image employed in example 1 of the method of the present invention;
FIG. 3 is a schematic diagram of a Raman spectrum employed in example 1 of the method of the present invention;
FIG. 4 is a schematic diagram of a smooth curve of the Raman spectrum obtained in example 1 of the method of the present invention;
FIG. 5 is a schematic diagram of the peak positions obtained in example 1 of the method of the present invention;
FIG. 6 is a schematic representation of a Gaussian fit image for parameter initialization according to method embodiment 1 of the invention;
FIG. 7 is a schematic representation of a Gaussian fit image after fine-tuning of the parameters of method embodiment 1 of the present invention;
FIG. 8 is a graph showing the Gaussian function of the 3 sigma range of method example 1 of the present invention.
FIG. 9 is a schematic diagram of the problem of unimodal overlap in example 1 of the process of the present invention.
FIG. 10 is a graphical representation of the Gaussian function after single peak de-duplication of example 1 of the method of the invention.
Detailed Description
The invention will be described in detail below with reference to the drawings and the detailed description.
Referring to fig. 1, the raman spectrum peak intensity calculation method of the present invention is implemented according to the following steps based on gaussian fitting principle:
step 1, acquiring an original spectrum image,
acquisition of spectral data in acquisition devicesReal number->Each element representing spectral data G, +.>Indicating that the spectral data G has m rows, n columns and d channels;
then, a spectral image is obtained after the operation according to the following formula (1)The expression is as follows:
where i denotes the ith row of the spectral data G and the spectral image H, j denotes the jth column of the spectral data G and the spectral image H, k denotes the kth channel of the spectral data G, and d is the total number of channels of the spectral data G.
It is apparent that the spectral image H has only one channel, and the spectral image H obtained in example 1 is shown in fig. 2.
Step 2, acquiring a Raman spectrum curve,
processing the spectrum image H obtained in the step 1 according to the following formula (2) to obtain a Raman spectrum curveThe Raman spectrum curve is the sum of the summation of each row of the spectrum image H,/and the like>Representation ofReal number,/->The expression for the raman spectral curve q, representing a vector of dimension n, is as follows:
the raman spectrum obtained in example 1 is shown in fig. 3.
Step 3, obtaining a Raman spectrum smooth curve,
for the raman spectrum curve obtained in step 2The raman spectrum smoothing curve q 'is obtained by processing according to the following formula (3), and the expression of the raman spectrum smoothing curve q' is as follows:
q′＝q*w (3)
wherein, the convolution is represented by,span represents the width of the smoothing window and the raman spectrum smoothing curve obtained in example 1 is shown in fig. 4.
Step 4, automatically obtaining the peak dividing positions,
4.1 Processing the Raman spectrum smoothing curve q' obtained in the step 3 according to the following formula (4) to obtain a mean valueMean->The expression of (2) is as follows:
the mean value obtained in example 1See the dashed line in fig. 4.
4.2 In the peak-splitting of the raman spectrum smoothing curve, a set Rs of the vertex x-axis values of each single peak is obtained,
for convenience of description, it is assumed that the peak of the raman spectrum smoothing curve includes L single peaks, that is, the peak x-axis values of the L single peaks included in the set Rs are denoted as Rs (L), l=1, 2, …, L, and the expression of the set Rs is calculated according to the following formula (5):
the peak-separated top x-axis values obtained in example 1 are shown in fig. 5 as black dots.
4.3 In the peak division of the raman spectrum smoothing curve, obtaining the x-axis pixel value range of a single peak,
step 4.2) obtaining values Rs (L) of the x-axis of the L unimodal vertices, l=1, 2, …, L, for the first unimodal, using values q (j) of the y-axis of the unimodal curve, j=1, 2, …, n being greater than the mean valueCan obtain the range of the first unimodal x-axis pixel value [ b (l), e (l) ]]，
The peak splitting ranges obtained in example 1 are shown in the graph of fig. 5.
Step 5, obtaining the physical value of the X-axis of the Raman spectrum curve,
since the values of the x-axis of the raman spectrum curve are only represented in pixel sense and do not contain physical sense, typically three sets of values are manually calibrated, corresponding to the pixel values and the physical values, so that the x-axis values of the raman spectrum curve have physical sense.
For convenience of description, three sets of pixel values and corresponding three sets of physical values are input, and the three sets of pixel values are defined asDefine the three sets of physical values as +.>For t pix And t phy Calculating according to the formula (6) to obtain a parameter matrix Q, wherein the expression of the parameter matrix Q is as follows:
based on parameter matrices Q and t phy The parameter vector c, c= [ c ] is calculated according to the following equation (7) 1 ,c 2 ,c 3 ]，c 1 ,c 2 ,c 3 Are all constants, the expression of the parameter vector c is as follows,
c＝t phy ·Q -1 (7)
the pixel value of the x-axis defining the raman spectrum curve is x pix ＝[1,2,…,n]The physical value x of the x-axis is calculated according to the following formula (8), and the expression of the physical value x is as follows:
step 6, initializing peak-splitting Gaussian parameters,
from the values Rs (L) of the x-axis of the vertices of the L unipeaks obtained in step 4.2), l=1, 2, …, L, and step 4.3), the x-axis ranges of the L unipeaks are obtained [ b (L), e (L)]L=1, 2, …, L, initializing gaussian parameters of L single peaks, including mean μ (L), standard deviation σ (L), maximum y max (l)，
6.1 Acquiring the discrete point coordinates of each single peak,
for the first single peak, according to its x-axis range [ b (l), e (l)]Obtaining an x-axis physical value vector x l ＝[x b(l) ,…,x e(l) ]Acquiring a y-axis physical value vector q l ＝[q b(l) ,…,q e(l) ]；
6.2 Acquiring the maximum value and the average value of each single peak,
initializing the mean μ (l) =x (Rs (l)) of the i-th unimodal x-axis physical value, wherein the value Rs (l) has been obtained in step 4.2); the x-axis physical value x (j), j=1, 2, …, n has been obtained in step 5;
initializing maximum value y max (l) Calculated according to equation (9), the expression is as follows:
let the discrete point of the first single peak have n l The number is:
6.3 A) an average residual value is obtained,
first, the standard deviation sigma is initialized 0 (l) =200, standard deviation σ 1 (l)＝10；
Then according to sigma 0 (l) Sum sigma 1 (l) The standard deviation, the average residual value of the Gaussian function value and the unimodal y-axis true value are calculated according to the following formulas (10) and (11), and the expressions are as follows:
wherein res 0 Is based on standard deviation sigma 0 (l) Calculating the Gaussian function value, wherein the Gaussian function value is an average residual value of the real value of the unimodal y-axis; res 1 Is based on standard deviation sigma 0 (l) Calculating the Gaussian function value, wherein the Gaussian function value is an average residual value of the real value of the unimodal y-axis; x is x s ∈[x b(l) ,…,x e(l) ]X-axis physical value expressed as the first single peak; q s ∈[q b(l) ,…,q e(l) ]Y-axis true value expressed as the first single peak; n is n l The number of discrete points expressed as the first single peak;
6.4 A) the average residual value is updated,
comparing the average residual value res obtained by the formulas (10) and (11) 0 And res 1 ，
If res 0 <res 1 Then update sigma 1 (l)＝0.2·σ 0 (l)+0.8·σ 1 (l) And update res using equation (11) 1 ；
Otherwise, update sigma 0 (l)＝0.8·σ 0 (l)+0.2·σ 1 (l) And update res using equation (10) 0 ；
6.5 A) each single peak standard deviation is obtained,
if res 0 >ε, and res 1 >Epsilon, go to step 6.4), where the parameter epsilon = 0.1;
otherwise, the standard deviation σ (l) = (σ) is obtained 0 (l)+σ 1 (l))/2。
The gaussian function curve after parameter initialization of example 1 is shown in fig. 6.
Step 7, fine tuning the split peak Gaussian parameters,
7.1 Calculating each single peak mean and corresponding residual mean,
first, mu is initialized around the mean mu (l) 0 (l) Sum mu 1 (l) Two means, mean μ 0 (l) Mu (l) -delta mu, mean mu 1 (l)＝μ(l)+Δμ；
Then, mu is utilized 0 (l) Sum mu 1 (l) Two average values are respectively obtained, and the average residual value of the two Gaussian function values and the unimodal y-axis true value is calculated to obtain a new average residual value diff 0 And diff (diff) 1 The expressions are respectively the following expression (12) and expression (13):
if diff 0 <diff 1 Then update mu 1 (l)＝0.2·μ 0 (l)+0.8·μ 1 (l) And update diff using equation (13) 1 ；
Otherwise, update mu 0 (l)＝0.8·μ 0 (l)+0.2·μ 1 (l) And update diff using equation (12) 0 ；
7.2 Fine-tuning each of the unimodal means,
if diff 0 >Epsilon' or diff 1 >Epsilon', go to step 7.1); wherein, parameter ε' =0.1;
otherwise, the mean μ (l) = (μ) is obtained 0 (l)+μ 1 (l))/2；
7.3 Fine-tuning the standard deviation of the standard,
initializing standard deviation sigma 0 (l) =2·σ (l), standard deviation σ 1 (l)＝σ(l)/2；
Then, the standard deviation σ (l) is updated according to step 6.5).
The gaussian function curve after fine tuning of the parameters of example 1 is shown in fig. 7.
Step 8, obtaining the peak splitting area of the Raman spectrum,
8.1 Initializing the unimodal x-axis range,
according to the 3σ principle, the probability of the x-axis range of the gaussian distribution being (μ -3σ, μ+3σ) is as high as 0.9974, and therefore, for each single peak, the x-axis rangeThe calculation is carried out according to the following formulas (14) and (15), and the formulas are as follows:
obtaining the x-axis range of each single peak according to the formulas (14) and (15)Thereafter, the maximum value y obtained according to step 6 max (l) And step 7 to obtainThe obtained mean μ (l) and standard deviation σ (l) are shown in fig. 8.
As can be seen from fig. 8, there is a phenomenon in which two peaks at the left end overlap. For clarity of illustration of the overlap problem, FIG. 9 is a visual depiction of the overlap problem for two single peaks.
To solve the problem of unimodal overlap as shown in FIG. 8 of example 1, i.eThus, a subsequent deduplication operation is also required, as follows:
8.2 A single peak x-axis range de-duplication,
initializing parametersWherein x is 0 <x 1 ，
According to parameter x 0 And x 1 The mean value x' is calculated according to the following formula (16) as follows:
x′＝(x 0 +x 1 )/2 (16)
obtaining Gaussian function value y by using Gaussian parameters and mean value x 0 Y 1 Calculated according to the following formulas (17) and (18), respectively, the expressions are as follows:
8.3 A) the parameters are updated and,
if it isUpdate->And turnTo step 8.4); wherein, parameter->
Otherwise, the Gaussian function value y obtained according to the step 8.2) is calculated 0 、y 1 Go through parameter x 0 、x 1 The update iteration operation of (a) is as follows:
if y 0 <y 1 Then update x 0 =x 'and update x' according to equation (16);
otherwise, update x 1 =x 'and update x' according to equation (16);
the gaussian function curve after single peak de-duplication of example 1 is shown in fig. 10.
8.4 The area of the peak of the raman spectrum is obtained,
according to the L single peaks obtained in the step 4.2), calculating L single peak areas S (L) according to the formula (19), wherein the expression of the single peak areas S (L) is as follows:
wherein,for the x-axis range of the first single peak, μ (l) is the mean parameter after the first single peak Gaussian fitting, σ (l) is the standard deviation parameter after the first single peak Gaussian fitting, y max (l) And (5) obtaining the maximum value parameter after the first single-peak Gaussian fitting.
The working principle of the invention is that firstly, the spectrum data collected by the equipment is converted into a Raman spectrum curve, and the number of the split peaks and the range of the single peak are automatically obtained in the Raman spectrum curve by utilizing a self-adaptive method; secondly, converting the x-axis pixel information of the Raman spectrum curve into physical information by using a polynomial fitting method; again, two steps of parameter initialization and parameter fine tuning are utilized to obtain the mean μ (l), standard deviation σ (l) and maximum y of different sub-peaks max (l) Three Gauss ginsengA number. Finally, the proposed iterative solution method is utilized to solve the problem of single-peak overlapping, and the areas of a plurality of single peaks are automatically obtained.
Example 1
The method according to the invention is carried out by the spectral image H obtained in step 1, see FIG. 2. The raman spectrum curve q obtained in step 2 is shown in fig. 3. The raman spectrum smoothing curve obtained in step 3 is shown in fig. 4. Step 4.1) the mean value of the raman spectrum smoothing curve is shown in fig. 4 by the dashed line. The x-axis value of the peak-splitting top point obtained in step 4.2) is shown in the black dot in fig. 5. The peak splitting ranges obtained in step 4.3) are shown in the graph of fig. 5. The gaussian function curve after the parameter initialization in step 6 is shown in fig. 6. The gaussian function curve after the fine adjustment of the parameters in step 7 is shown in fig. 7. Step 8.1) initializing a gaussian function curve for the unimodal x-axis range as shown in fig. 8. The single peak overlap problem present in fig. 8 is visually depicted in fig. 9. Step 8.2) the gaussian function curve after unimodal de-duplication is shown in fig. 10.
The implementation parameters and results are: the implementation process refers to the method of example 1, and the method of automatically obtaining the peak-dividing peak number in step 4 lays a foundation for automatically obtaining peak-dividing area in the invention; step 5 requires manual input of three sets of pixel valuesAnd corresponding three sets of physical values +.>The invention automatically acquires the physical value of the X-axis of the Raman spectrum curve. Step 6 and step 7 finally obtain three gaussian parameters of different sub-peaks: mean μ (k), standard deviation σ (k), maximum y max (k) A. The invention relates to a method for producing a fibre-reinforced plastic composite Step 8 is to solve the problem of unimodal overlap and automatically obtain the areas of multiple unipeaks.
Table 1 shows the experimental results of example 1, with rows 1 "Mono 1", "Mono 2", "Mono 3", "Mono 4" corresponding to four Mono peaks in the order from left to right in FIG. 10. The calculated area value of the second row is the calculation result of the algorithm of the invention; the third row is the actual area calculated by the value of the actual curve and the unimodal x-axis physical value range obtained in the step 8.3); the accuracy of the four single peaks is 68.83 percent at the minimum; up to 91.42% and an average accuracy of up to 82.75%. In summary, the invention can automatically calculate the peak splitting intensity of the Raman spectrum without manually determining the x-axis range of the single peak.
Table 1, example 1 Experimental results
| Single peak 1 | Single peak 2 | Single peak 3 | Single peak 4 | |
| Calculation area (10) 8 *) | 0.2807 | 1.6650 | 1.0916 | 0.9609 |
| Actual area (10) 8 *) | 0.4078 | 1.8212 | 1.2707 | 1.1327 |
| Accuracy (%) | 68.83 | 91.42 | 85.91 | 84.83 |
Example 2
Implementation procedure referring to example 1, the implementation parameters and results are: step 6 of the invention firstly carries out the mean value mu (l), the standard deviation sigma (l) and the maximum value y of three Gaussian parameters max (l) And (3) obtaining a unimodal maximum value and a unimodal average value in an effective range, and obtaining a standard deviation in an iterative mode. The calculation complexity of step 6 is O (n), and the least squares method is O (n 2 ) The method of the invention is superior because of the large amount; step 7 of the invention is carried out to three Gaussian parameter mean values mu (l), standard deviation sigma (l) and maximum value y max (l) Fine tuning of (c) to obtain more accurate parameters.
As shown in table 2, the experimental comparison results of the deviation between the calculated function value and the true value are shown in the table, and the three data in the table are the related experimental results of the least square method, the method of step 6 of the present invention and the method of the present invention (step 6+step 7), respectively. Since the raman spectrum peak separation is only similar to a gaussian function, an error may exist in the acquired single peak range, and the gaussian function parameter is obtained by using a least square method according to the data of the range, and the deviation between the function value and the true value is shown in table 2, and the deviation result of the least square method is far higher than the deviation of the invention. As can be seen from the experimental data of the four single peaks, the deviation of the present invention is smaller than that of step 6. Therefore, the Gaussian parameter fine adjustment method of the step 7 can improve the parameter precision.
Table 2, experimental results of the deviation of the function value from the true value in example 2
| Single peak 1 | Single peak 2 | Single peak 3 | Single peak 4 | |
| Least squares method (10 5 *) | 0.4391 | 6.7227 | 5.3613 | 6.1244 |
| Step 6 method (10) 5 *) | 0.3155 | 1.3274 | 1.4427 | 1.4091 |
| The invention (10) 5 *) | 0.2562 | 1.3048 | 1.3038 | 1.3266 |
Example 3
Implementation procedure referring to example 1, the implementation parameters and results are: step 8 of the present invention automatically obtains the areas of the plurality of single peaks. First, according to the 3 sigma principle, the x-axis range of the single peak is automatically obtained, as shown in table 3, and the overlapping phenomenon occurs between the right end of the single peak 1 and the left end of the single peak 2. According to the de-overlapping method of step 8.3) of the present invention, it is calculated that the right end of the single peak 1 and the left end of the single peak 2 are 1319.2. The invention solves the problem of single peak overlapping and lays a foundation for accurately obtaining the single peak area subsequently.
TABLE 3 automatic acquisition of unimodal x-axis Range results from step 8
| Single peak 1 | Single peak 2 | Single peak 3 | Single peak 4 | |
| Left end (10) 3 *) | 1.0297 | 1.2824 | 2.6458 | 3.4344 |
| Right end (10) 3 *) | 1.3876 | 1.5623 | 2.9148 | 3.6738 |
In conclusion, the method can solve the defects in the prior art and realize accurate calculation of the peak splitting intensity of the Raman spectrum.
Claims (9)
1. The Raman spectrum peak splitting intensity calculating method based on Gaussian fitting is characterized by comprising the following steps of:
step 1, acquiring an original spectrum image;
step 2, acquiring a Raman spectrum curve;
step 3, acquiring a Raman spectrum smooth curve;
step 4, automatically obtaining peak dividing positions;
step 5, acquiring an x-axis physical value of a Raman spectrum curve;
step 6, initializing peak-splitting Gaussian parameters;
step 7, fine tuning the split peak Gaussian parameters;
and 8, acquiring the peak splitting area of the Raman spectrum.
2. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 1, the specific process is as follows:
acquiring spectral dataReal number->Each element representing spectral data G, +.>Indicating that the spectral data G has m rows, n columns and d channels;
then, a spectral image is obtained after the operation according to the following formula (1)The expression is as follows:
where i denotes the ith row of the spectral data G and the spectral image H, j denotes the jth column of the spectral data G and the spectral image H, k denotes the kth channel of the spectral data G, and d is the total number of channels of the spectral data G.
3. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 2, the specific process is as follows:
processing the spectrum image H obtained in the step 1 according to the following formula (2) to obtain a Raman spectrum curveThe Raman spectrum curve is the sum of the summation of each row of the spectrum image H,/and the like>Representing real number,/->The expression for the raman spectral curve q, representing a vector of dimension n, is as follows:
4. the method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 3, the specific process is as follows:
for the raman spectrum curve obtained in step 2The raman spectrum smoothing curve q 'is obtained by processing according to the following formula (3), and the expression of the raman spectrum smoothing curve q' is as follows:
q′＝q*w (3)
wherein, the convolution is represented by,span represents the width of the smoothing window.
5. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 4, the specific process is as follows:
4.1 Processing the Raman spectrum smoothing curve q' obtained in the step 3 according to the following formula (4) to obtain a mean valueMean->The expression of (2) is as follows:
4.2 In the peak-splitting of the raman spectrum smoothing curve, a set Rs of the vertex x-axis values of each single peak is obtained,
assuming that the peak of the raman spectrum smoothing curve contains L single peaks, i.e. the peak x-axis values of the L single peaks contained in the set Rs are denoted as Rs (L), l=1, 2, …, L, and the expression of the set Rs is calculated according to the following formula (5):
4.3 In the peak division of the raman spectrum smoothing curve, obtaining the x-axis pixel value range of a single peak,
step 4.2) obtains the values Rs (L) of the x-axis of the L unimodal vertices, l=1, 2, …, L, and for the L-th unimodal, the range of pixel values [ b (L), e (L) ] of the L-th unimodal vertex can be obtained using the constraint that the value q (j), j=1, 2, …, n of the y-axis of the unimodal curve is greater than the mean value q'.
6. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 5, the specific process is as follows:
three groups of pixel values and corresponding three groups of physical values are input, and the three groups of pixel values are defined asDefine the three sets of physical values as +.>For t pix And t phy Calculating according to the formula (6) to obtain a parameter matrix Q, wherein the expression of the parameter matrix Q is as follows:
based on parameter matrices Q and t phy The parameter vector c, c= [ c ] is calculated according to the following equation (7) 1 ,c 2 ,c 3 ]，c 1 ,c 2 ,c 3 Are all constants, the expression of the parameter vector c is as follows,
c＝t phy ·Q -1 (7)
the pixel value of the x-axis defining the raman spectrum curve is x pix ＝[1,2,…,n]The physical value x of the x-axis is calculated according to the following formula (8), and the expression of the physical value x is as follows:
7. the method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 6, the specific process is as follows:
from the values Rs (L) of the x-axis of the vertices of the L unipeaks obtained in step 4.2), l=1, 2, …, L, and step 4.3), the x-axis ranges of the L unipeaks are obtained [ b (L), e (L)]L=1, 2, …, L, initializing gaussian parameters of L single peaks, including mean μ (L), standard deviation σ (L), maximum y max (l)，
6.1 Acquiring the discrete point coordinates of each single peak,
for the first single peak, according to its x-axis range [ b (l), e (l)]Obtaining an x-axis physical value vector x l ＝[x b(l) ,…,x e(l) ]Acquiring a y-axis physical value vector q l ＝[q b(l) ,…,q e(l) ]；
6.2 Acquiring the maximum value and the average value of each single peak,
initializing the mean μ (l) =x (Rs (l)) of the i-th unimodal x-axis physical value, wherein the value Rs (l) has been obtained in step 4.2); the x-axis physical value x (j), j=1, 2, …, n has been obtained in step 5;
initializing maximum value y max (l) Calculated according to equation (9), the expression is as follows:
let the discrete point of the first single peak have n l The number is:
6.3 A) an average residual value is obtained,
first, the standard deviation sigma is initialized 0 (l) =200, standard deviation σ 1 (l)＝10；
Then according to sigma 0 (l) Sum sigma 1 (l) The standard deviation, the average residual value of the Gaussian function value and the unimodal y-axis true value are calculated according to the following formulas (10) and (11), and the expressions are as follows:
wherein res 0 Is based on standard deviation sigma 0 (l) Calculating the Gaussian function value, wherein the Gaussian function value is an average residual value of the real value of the unimodal y-axis; res 1 Is based on standard deviation sigma 0 (l) Calculated Gaussian functionA value, which is an average residual value of the function value and the unimodal y-axis true value; x is x s ∈[x b(l) ,…,x e(l) ]X-axis physical value expressed as the first single peak; q s ∈[q b(l) ,…,q e(l) ]Y-axis true value expressed as the first single peak; n is n l The number of discrete points expressed as the first single peak;
6.4 A) the average residual value is updated,
comparing the average residual value res obtained by the formulas (10) and (11) 0 And res 1 ，
If res 0 <res 1 Then update sigma 1 (l)＝0.2·σ 0 (l)+0.8·σ 1 (l) And update res using equation (11) 1 ；
Otherwise, update sigma 0 (l)＝0.8·σ 0 (l)+0.2·σ 1 (l) And update res using equation (10) 0 ；
6.5 A) each single peak standard deviation is obtained,
if res 0 >ε, and res 1 >Epsilon, go to step 6.4), where the parameter epsilon = 0.1;
otherwise, the standard deviation σ (l) = (σ) is obtained 0 (l)+σ 1 (l))/2。
8. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 7, the specific process is as follows:
7.1 Calculating each single peak mean and corresponding residual mean,
first, mu is initialized around the mean mu (l) 0 (l) Sum mu 1 (l) Two means, mean μ 0 (l) Mu (l) -delta mu, mean mu 1 (l)＝μ(l)+Δμ；
Then, mu is utilized 0 (l) Sum mu 1 (l) Two average values are respectively obtained, and the average residual value of the two Gaussian function values and the unimodal y-axis true value is calculated to obtain a new average residual value diff 0 And diff (diff) 1 The expressions are respectively the following expression (12) and expression (13):
if diff 0 <diff 1 Then update mu 1 (l)＝0.2·μ 0 (l)+0.8·μ 1 (l) And update diff using equation (13) 1 ；
Otherwise, update mu 0 (l)＝0.8·μ 0 (l)+0.2·μ 1 (l) And update diff using equation (12) 0 ；
7.2 Fine-tuning each of the unimodal means,
if diff 0 >Epsilon' or diff 1 >Epsilon', go to step 7.1); wherein, parameter ε' =0.1;
otherwise, the mean μ (l) = (μ) is obtained 0 (l)+μ 1 (l))/2；
7.3 Fine-tuning the standard deviation of the standard,
initializing standard deviation sigma 0 (l) =2·σ (l), standard deviation σ 1 (l)＝σ(l)/2；
Then, the standard deviation σ (l) is updated according to step 6.5).
9. The method for calculating the peak intensity of the raman spectrum based on gaussian fitting according to claim 1, wherein in step 8, the specific process is as follows:
8.1 Initializing the unimodal x-axis range,
x-axis range for each single peakThe calculation is carried out according to the following formulas (14) and (15), and the formulas are as follows:
obtaining the x-axis range of each single peak according to the formulas (14) and (15)Thereafter, the maximum value y obtained according to step 6 max (l) And the mean value mu (l) and the standard deviation sigma (l) obtained in the step 7;
to solve the problem of unimodal overlap, i.eSubsequent de-duplication operations are also required, as follows:
8.2 A single peak x-axis range de-duplication,
initializing parametersWherein x is 0 <x 1 ，
According to parameter x 0 And x 1 The mean value x' is calculated according to the following formula (16) as follows:
x′＝(x 0 +x 1 )/2 (16)
obtaining Gaussian function value y by using Gaussian parameters and mean value x 0 Y 1 Calculated according to the following formulas (17) and (18), respectively, the expressions are as follows:
8.3 A) the parameters are updated and,
if it isUpdate->And go to step 8.4); wherein, parameter->
Otherwise, the Gaussian function value y obtained according to the step 8.2) is calculated 0 、y 1 Go through parameter x 0 、x 1 The update iteration operation of (a) is as follows:
if y 0 <y 1 Then update x 0 =x 'and update x' according to equation (16);
otherwise, update x 1 =x 'and update x' according to equation (16);
8.4 The area of the peak of the raman spectrum is obtained,
according to the L single peaks obtained in the step 4.2), calculating L single peak areas S (L) according to the formula (19), wherein the expression of the single peak areas S (L) is as follows:
wherein,for the x-axis range of the first single peak, μ (l) is the mean parameter after the first single peak Gaussian fitting, σ (l) is the standard deviation parameter after the first single peak Gaussian fitting, y max (l) The maximum value parameter after the first single-peak Gaussian fitting is obtained.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311353099.9A CN117574064A (en) | 2023-10-18 | 2023-10-18 | Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311353099.9A CN117574064A (en) | 2023-10-18 | 2023-10-18 | Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN117574064A true CN117574064A (en) | 2024-02-20 |
Family
ID=89892511
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311353099.9A Pending CN117574064A (en) | 2023-10-18 | 2023-10-18 | Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117574064A (en) |
-
2023
- 2023-10-18 CN CN202311353099.9A patent/CN117574064A/en active Pending
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109800779B (en) | Change detection method for fusing FCM algorithm by using D-S evidence theory | |
| CN111709901B (en) | FCM cluster matching and Wallis filtering-based no-weight multi/hyperspectral remote sensing image color homogenizing method | |
| CN113762208B (en) | Spectrum conversion method of near infrared spectrum and characteristic spectrum and application thereof | |
| CN113177548B (en) | Key area identification method for immune fixed electrophoresis | |
| CN114323536B (en) | Interpolation method for improving measurement accuracy of five-hole probe | |
| CN116168037B (en) | Method and system for calculating bending degree of wire crimping based on image processing | |
| CN104598922B (en) | Full polarimetric SAR sorting technique based on fuzzy C-mean algorithm | |
| CN108921170B (en) | Effective image noise detection and denoising method and system | |
| CN112415078A (en) | Mass spectrum data spectrogram signal calibration method and device | |
| CN117574064A (en) | Gaussian fitting-based Raman spectrum peak-splitting intensity calculation method | |
| CN109300169B (en) | Semitransparent image color migration method based on linear transformation | |
| CN113359435B (en) | Correction method for dynamic working condition data of thermal power generating unit | |
| CN112329733B (en) | Winter wheat growth monitoring and analyzing method based on GEE cloud platform | |
| CN107578445B (en) | Image discriminable region extraction method based on convolution characteristic spectrum | |
| CN112052057B (en) | Data visualization method and system for optimizing color chart based on spring model | |
| CN113870130A (en) | Low-rank tensor completion method based on three-dimensional total variation and Tucker decomposition | |
| CN110705132B (en) | Guidance control system performance fusion evaluation method based on multi-source heterogeneous data | |
| CN107014785A (en) | A kind of improved method of emission spectrum background correction | |
| CN107705244B (en) | Edge connection correction method suitable for large-area multi-remote sensing image | |
| CN112666097B (en) | Fourier infrared polynomial iterative smoothing filtering baseline correction method | |
| CN109919872A (en) | A kind of image recovery method, system, readable storage medium storing program for executing and computer equipment | |
| WO2023201772A1 (en) | Cross-domain remote sensing image semantic segmentation method based on adaptation and self-training in iteration domain | |
| CN113223098B (en) | Preprocessing optimization method for image color classification | |
| CN109903349B (en) | Color harmony prediction method based on maximum likelihood estimation | |
| CN111126452A (en) | Ground feature spectral curve expansion method and system based on principal component analysis |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |