An arbitrary-order calculus method for spectra

By separating signal and noise in the frequency domain using the SH operator and implementing arbitrary-order calculus in the time domain using Fourier transform and convolution kernel, the problem of low noise filtering rate and signal loss of differentiators in noisy environments in existing technologies is solved, and high-precision chemical spectral signal processing is achieved.

CN116306960BActive Publication Date: 2025-11-14SUZHOU DALAI SOFTWARE TECH CO LTD
View PDF 2 Cites 0 Cited by

Patent Information

Application Number
CN202211685779.6
Authority / Receiving Office
CN · China
Patent Type
Patents(China)
Current Assignee / Owner
Filing Date
2022-12-27
Publication Date
2025-11-14
Estimated Expiration
2042-12-27

AI Technical Summary

Technical Problem

Existing technologies lack high-order calculus tools suitable for chemical analysis, especially in the presence of noise, which leads to low noise filtering efficiency of differentiators, signal loss, and problems with the confidence of results.

Method used

Arbitrary-order calculus is performed using the SH operator. By simulating hyperbolic functions, signals and noise are separated in the frequency domain. Arbitrary-order calculus of signals is realized in the time domain using Fourier transform and convolution kernels. The SH operator (SHfd) is designed to be executed in MATLAB.

Benefits of technology

It achieves accurate signal separation and high-precision arbitrary-order calculus in noisy environments, improves the robustness and accuracy of chemical spectral signal processing, and is suitable for calculating integer and fractional derivatives.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN116306960B_ABST
    Figure CN116306960B_ABST
Patent Text Reader

Abstract

An arbitrary-order calculus method for spectra is proposed. Based on the properties of the data, a more accurate calculus operator applicable to arbitrary orders is designed by simulating hyperbolic functions. This method can perform arbitrary-order calculus on various signals, including integer and fractional derivatives and integrals. The method is computationally simple and possesses robustness and the necessary accuracy.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] An arbitrary-order calculus method for spectra Background Technology

[0002] In analytical chemistry, derivative spectroscopy (DS) is fundamental for handling overlapping spectra and baseline correction, and is a commonly used method to improve selectivity, in addition to multivariate statistical methods. Despite decades of development, DS still has some unresolved issues, such as the lack of suitable calculus calculators.

[0003] Noise-robust differentiators are rare, even under very weak noise conditions. Measurements cannot always be free of noise, and there has been a lack of chemist-friendly, practical calculus tools, especially higher-order differentiators, suitable for real-world data. Before the advent of higher-order derivatives, almost all solvers were intolerant of noise and lacked accuracy.

[0004] Reported methods include wavelets, Fourier transforms, and finite-difference gradient denoising, with the SG method being the most widely accepted. Originating in 1964, the SG differentiator combines polynomial fitting and linear least squares. It has been widely accepted by analytical chemists due to its balance between denoising efficiency and ease of use, particularly its acceptability at lower derivatives. Even so, at higher orders, the SG differentiator suffers from low noise filtering efficiency, signal loss, and issues with the confidence level of the results. A user-friendly calculator for chemists should be powerful, accurate, and easy to use. Currently, derivative spectroscopy still suffers from limitations due to the lack of suitable tools. Summary of the Invention

[0005] To address the problems existing in the prior art, the present invention aims to provide an arbitrary-order calculus method for spectra;

[0006] To achieve the above objectives, the technical solution of the present invention is as follows:

[0007] An arbitrary-order calculus method for spectra consists of the following steps:

[0008] 1. Spectral acquisition to obtain the raw spectrum to be processed.

[0009] The sample spectrum can be collected, but is not limited to visible, ultraviolet, fluorescence, Raman, near-infrared, mid-infrared, etc.

[0010] 2. Using the SH operator, calculate the arbitrary-order calculus of the spectrum.

[0011] The measured spectral signal contains both signal and noise, i.e., the measured value = signal + noise. In other words, in the frequency domain, the envelope of the measured value is the sum of the signal envelope and the noise envelope.

[0012] M = S + N(1)

[0013] For white noise, N is a constant n, so S is the product of the known M and the response function K.

[0014] S = M·K(2)

[0015] So

[0016]

[0017] Since the real number field is discontinuous, it can be simplified by simulating a hyperbolic function.

[0018] Redefining K as

[0019] SH=asech(bω) p (4)

[0020] Formula (4) includes three parameters: amplitude a, width (for stopping frequency) b, and peak steepness exponent p. This function combines the constant n into p, which also means that even for the same signal, the response function corresponding to different levels of noise is different, and the numerical setting can be omitted.

[0021] Equation (5) separates the signal in the frequency domain and then performs an inverse Fourier transform to transfer it to the time domain.

[0022]

[0023] Using convolution for time-domain measurements

[0024]

[0025]

[0026] k is the noise reduction convolution kernel.

[0027] So

[0028] s=m*k (8)

[0029] According to the properties of convolution, the νth derivative of the signal is:

[0030] s (v) =m*k (v) (9)

[0031] This property applies to real numbers ν. In other words, any order of calculus can be implemented on a convolution kernel with sufficient accuracy.

[0032] The SH function mentioned earlier has three parameters, of which the amplitude parameter 'a' is not mandatory. To maintain the same magnitude before and after convolution, the kernel k has a necessary condition:

[0033]

[0034] Therefore, only b and p are necessary in the Sech function.

[0035] invSH(b,p)=real(iFFT(sech(bω) p ))(11)

[0036] k = invSH / sum(invSH)(12)

[0037] k can undertake various tasks of calculus of any order.

[0038] A SH operator (S = SHfd(Mraw,ν,b,p,ON)) capable of performing arbitrary derivative calculus calculations was designed. It is an executable script in MATLAB. The sub-function Sechpf generates a time-domain convolution kernel based on the simulated frequency domain response Sech and the sub-function GLfd, which is used for GL fractional-order calculations.

[0039] 2) Input the original spectrum and set the parameters.

[0040] The calculation S = SHfd(Mraw, ν, b, p, ON) has 5 input terms, of which Mraw and the derivative order ν are required.

[0041] a) The first parameter, Mraw, is the data sequence to be processed, and the second parameter is the calculus order ν. When only the raw data Mraw and the derivative order ν are input, by default there is no signal filtering or denoising; it is simply a calculator for arbitrary derivatives. The order is not necessarily a positive integer; the derivative is positive, the integral is negative, and the fractional order is a non-integer.

[0042] b) The third input parameter is the peak width b, which is mainly adjusted based on the signal noise intensity.

[0043] c) The setting of the fourth parameter p determines the filter's adaptability to various frequencies, similar to the polynomial order in SG. This parameter does not need to be an integer, nor is there a restriction on the width setting. Fine-tuning this parameter can better approximate the true value, an improvement that greatly increases the flexibility of the method.

[0044] d) Any value can be entered in the fifth parameter, the purpose of which is to activate the drawing interface used for tuning.

[0045] 3) Run the program with the parameters set, and output the processed result S.

[0046] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0047] Chemical spectral signals tend to be low-frequency, while noise is uniformly distributed across the frequency spectrum. Whether integer or fractional, derivative operations can suppress low frequencies, but exhibit exponential enhancement at their respective orders. This contradiction leads to the derivative robustness problem of noisy signals. To address this issue, a more accurate representation of the spectrum in the data is needed. The three typical spectral peak modes—Lorentz, Voigt, and Gaussian—have consistent exponential function expressions in the frequency domain. When the multi-peak transformation from the time domain to the frequency domain is expressed as the sum of multiple exponential distributions, forming an envelope that preserves the overall distribution of the exponential expression. Therefore, an accurate signal response exists in the mixed spectrum of noise and signal. This invention proposes a hyperbolic sine function that can well simulate spectral signals.

[0048] By performing an inverse Fourier transform of the simulated hyperbolic function back to the time domain, a convolution kernel is obtained, which can then separate the signal from noisy data. Operations around the convolution kernel can then perform arbitrary-order calculus on the original spectrum.

[0049] Based on the properties of the data, this invention designs a more accurate arbitrarily ordered calculus operator (SH) for hyperbolic functions. Compared with existing estimation methods such as integer polynomial fitting, this operator can more accurately describe and approximate the actual distribution of the spectrum using hyperbolic response functions, making numerical calculations more accurate. Results from rigorously analytic functions verify the accuracy of the SH algorithm and its robustness to noisy signal processing.

[0050] The SH algorithm proposed in this invention can perform arbitrary-order differentials and integrals for various signals, including integer and fractional derivatives. This method is computationally very simple and possesses robustness and necessary accuracy. The applications of differentiation in analytical chemistry are countless; fractional derivative spectroscopy is helpful for multi-component quantification and peak fitting, among other things. Attached Figure Description

[0051] Figure 1 Fractional derivative of the sine function

[0052] Figure 2 Comparison of the accuracy results of the SH algorithm and the GL algorithm with noisy sine function

[0053] Figure 3 Optimize the interface diagram, among which Figure 3 (a)S=SHfd(LYS',2,40,2.6,0), Figure 3 (b)S=SHfd(LYS',2,25,1.55,0);

[0054] Figure 4 The original Raman spectrum of acetaminophen

[0055] Figure 5 Derivatives of integer orders (1-6) of the Raman spectrum of acetaminophen

[0056] Figure 6 Fractional derivatives of acetaminophen Raman spectra

[0057] Figure 7 Results of different orders of Raman spectroscopy for acetaminophen

[0058] Figure 8 Near-infrared spectra of sucrose, fructose, glucose, and raw brown sugar

[0059] Figure 9 Normalized derivative spectra of fructose from order 1 to 4

[0060] Figure 10 The fractional order of the glucose spectral data was obtained. Figure 10 (a) shows the results of normalizing the fractional derivatives for all data. Figure 10 (b) is a magnified view of a portion of the image.

[0061] Figure 11 Results of 1.8th order integrals for sucrose, fructose, glucose, and raw brown sugar. Detailed Implementation

[0062] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments: For example... Figure 1-11 As shown,

[0063] The SH method is a numerical computation tool with advantages such as robustness, accuracy, ease of use, and no need for cumbersome setup. While fractional calculus has many of its own properties, it is, in principle, calculus and must satisfy the fundamental definition of calculus itself. Like integer calculus, it is a linear operator, satisfying all the basic arithmetic rules of everyday chemical calculations, such as the commutative, associative, and distributional laws. Therefore, chemists do not need to get bogged down in lengthy and tedious specific definitions. In fact, calculus is ubiquitous in chemical data processing and preprocessing, not only in derivatives and spectroscopy.

[0064] An arbitrary-order calculus method for spectra consists of the following steps:

[0065] 1. Spectral acquisition to obtain the raw spectrum to be processed.

[0066] The sample spectrum can be collected, but is not limited to visible, ultraviolet, fluorescence, Raman, near-infrared, mid-infrared, etc.

[0067] 2. Using the SH operator, calculate the arbitrary-order calculus of the spectrum.

[0068] The measured spectral signal contains both signal and noise, i.e., the measured value = signal + noise. In other words, in the frequency domain, the envelope of the measured value is the sum of the signal envelope and the noise envelope.

[0069] M = S + N (1)

[0070] For white noise, N is a constant n, so S is the product of the known M and the response function K.

[0071] S = M·K (2)

[0072] So

[0073]

[0074] Since the real number field is discontinuous, it can be simplified by simulating a hyperbolic function.

[0075] Redefining K as

[0076] SH=asech(bω) p (4)

[0077] Formula (4) includes three parameters: amplitude a, width (for stopping frequency) b, and peak steepness exponent p. This function combines the constant n into p, which also means that even for the same signal, the response function corresponding to different levels of noise is different, and the numerical setting can be omitted.

[0078] Equation (5) separates the signal in the frequency domain and then performs an inverse Fourier transform to transfer it to the time domain.

[0079]

[0080] Using convolution for temporal measurements

[0081]

[0082]

[0083] k is the noise reduction convolution kernel.

[0084] So

[0085] s=m*k (8)

[0086] According to the properties of convolution, the νth derivative of the signal is:

[0087] s (v) =m*k (v) (9)

[0088] This property applies to real numbers ν. In other words, any order of calculus can be implemented on a convolution kernel with sufficient accuracy.

[0089] The SH function mentioned earlier has three parameters, of which the amplitude parameter 'a' is not mandatory. To maintain the same magnitude before and after convolution, the kernel k has a necessary condition:

[0090]

[0091] Therefore, only b and p are necessary in the Sech function.

[0092] invSH(b,p)=real(iFFT(sech(bω) p (11)

[0093] k = invSH / sum(invSH) (12)

[0094] k can undertake various tasks of calculus of any order.

[0095] A SH operator (S = SHfd(Mraw,ν,b,p,ON)) capable of performing arbitrary derivative calculus calculations was designed. It is an executable script in MATLAB. The sub-function Sechpf generates a time-domain convolution kernel based on the simulated frequency domain response Sech and the sub-function GLfd, which is used for GL fractional-order calculations.

[0096] 2) Input the original spectrum and set the parameters.

[0097] The computation S = SHfd(Mraw, ν, b, p, ON) has 5 input terms, of which Mraw and the derivative order ν are required.

[0098] a) The first parameter, Mraw, is the data sequence to be processed, and the second parameter is the calculus order ν. When only the raw data Mraw and the derivative order ν are input, by default there is no signal filtering or denoising; it is simply a calculator for arbitrary derivatives. The order is not necessarily a positive integer; the derivative is positive, the integral is negative, and the fractional order is a non-integer.

[0099] b) The third input parameter is the peak width b, which is mainly adjusted based on the signal-noise intensity.

[0100] c) The setting of the fourth parameter p determines the filter's adaptability to various frequencies, similar to the polynomial order in SG. This parameter does not need to be an integer, nor is there a restriction on the width setting. Fine-tuning this parameter can better approximate the true value, an improvement that greatly increases the flexibility of the method.

[0101] d) Any value can be entered in the fifth parameter, the purpose of which is to activate the drawing interface used for tuning.

[0102] 3) Run the program with the parameters set, and output the processed result S.

[0103] Example 1

[0104] 1) Discrete points are selected within the step size h. The normalized Gaussian peak is then analyzed using MATLAB's built-in diff and Gradient functions, the minimum window SG method, the GL numerical algorithm, and the SH method proposed in this invention. Perform 2 nd -10 th Find the derivative.

[0105] 2) The SH parameters are set to b=35, p=7.8, and the SG parameters are set to width=11, order=2. Other algorithms are directly called. The phase shift Δstep in the X direction and the truncation error in the Y direction are used to calculate the results. The accuracy of the numerical derivation of the five methods was examined, and the results are shown in Table 1.

[0106] Table 1. Comparison of second- to tenth-order accuracy between the SH algorithm and commonly used derivative algorithms.

[0107]

[0108]

[0109] **Mathematical analytical formulas for derivatives of all orders, where h is the step size.

[0110] The results show that the phase shifts of diff and GL are forward and backward shifts, respectively, and increase with increasing order. SH, Gradient, and SG can maintain no drift. After the 8th derivative, except that the truncation error of diff is 4 times larger than that of the original method, the truncation error of the SH method is consistent with that of other methods.

[0111] Due to floating-point errors, neither Diff nor GL can produce stable results when the order is increased to the ultra-high order of 10. Considering that the error is extremely weak noise, SG and SH can obtain stable results by further increasing the window width on the basis of the 8th order.

[0112] The results show that the SH method achieves the highest accuracy across different derivatives and has passed rigorous mathematical analysis. It exhibits the lowest truncation error without phase shift. Even under extreme conditions, it meets the requirements of robustness, accuracy, and convenience.

[0113] Example 2

[0114] 1) Sin function

[0115] The calculation of the fractional derivative includes the phase shift in the X direction and its value in the Y direction. The sine function has a quantitative relationship with the order and magnitude of the phase shift, and is the most commonly used evaluation criterion for algorithms. The fractional derivative of the sine function is as follows: Figure 1 As shown.

[0116] The derivative property of Sin(x) is

[0117]

[0118] h is the step size, and v is the order. We have:

[0119]

[0120]

[0121] For the coordinate tracking of a peak, equations (2) and (3) can provide accurate differential results. 2) Comparison of the results of SH and GL numerical methods for processing the sine function

[0122] (1) Comparison of processing results under noise-free conditions

[0123] When v < 8, the GL method maintains the linear relationship between X and v. However, phase drift exhibits cumulative error, with the 8th order phase drift exceeding 4 steps. The regression equation for X and ν is X = / 2.03ν, RMSE = 5.0e -3 For regressions below the 8th order, h = pi / 100, the regression equation is Ln(Y) / Ln(h) = ν, and RMSE = 1.36. -5 .

[0124] Without setting parameters b and p, SH maintains a good linear relationship between X and ν until ν > 8, with the regression equation being X = pi / 2ν and RMSE = 5.1e. -14 The regression equation for Y with h = pi / 100 is Ln(Y) / Ln(h) = ν, RMSE = 2.5e -6 The tiny error comes from floating-point numbers. The calculation results of SH are in perfect agreement with the analytical equations, indicating that its accuracy in both X and Y is better than that of GL. The floating-point error is only half that of GL. More importantly, SH does not have phase drift, which is beneficial for maintaining an accurate description of the spectral wavelength.

[0125] (2) Comparison of processing results under noisy conditions

[0126] Noise will compromise the robustness of differentials, whether they are integer or fractional. Figure 2 A comparison of the computational accuracy results of the SH algorithm and the GL algorithm with a noisy sine function. Figure 2 It can be seen that after adding one-thousandth of white noise, the GL algorithm can no longer maintain the linear relationship between X and ν after the second derivative.

[0127] When SH(original data Mraw, ν, 350, 2) is set, the results of the SH algorithm at any order up to the 8th order still strictly conform to the formula solution, with RMSEs of 5.1e for X and Y respectively. -14 and 3.5e-4 .

[0128] Example 3

[0129] 1) Collect the ultraviolet spectrum of chlorogenic acid (UV-Vis spectrophotometer, S3150, Scinco, South Korea), i.e., raw data LYS.

[0130] 2) In S = SHfd(Mraw, v, b, p, ON), input the original spectral data and corresponding parameters, run the program, and output the tuning interface as shown below. Figure 3 , Figure 3 (a) Figure 3 (b) shows the calculation results for different parameters.

[0131] The top left corner shows a comparison between the original data and the data after 0th-order processing, while the top right corner shows the output after 2nd-order differentiation. The two graphs below show the residuals before and after differentiation, i.e., the filtered noise; the left graph represents the time domain, and the right graph represents the frequency domain. Visualizing the results can help adjust the settings of parameters (ν, b, p). If the noise in the two graphs is relatively flat, it indicates that the filtering is appropriate. In the frequency graph, if the middle part shrinks, it indicates that the low-frequency residuals are significant; if there are bulges, it indicates that over-denoising has occurred.

[0132] Example 3

[0133] 1) A handheld Raman spectrometer (Platinum Jingtan NX10 spectrometer, sampling pixel 2048, spectral resolution 8cm²) was used. -1 Spectral range 300-3000cm -1 At an integration time of 100 ms, the Raman spectrum of acetaminophen (AR) was acquired to obtain the raw spectrum Mraw. Data points from 180 to 1000 were extracted to obtain X, as shown in the spectrum. Figure 4 As shown.

[0134] 2) Set SH(original data X, νv, 15, 2.6, 0), where v is of order 1-6, with an interval of 1 order. See the results. Figure 5 .

[0135] 3) Set SH(original data X, ν, 15, 2.6), where v is of order 1.1-2 with an interval of order 0.1. The result after normalizing the derivative spectrum is shown in [the original data X, ν, 15, 2.6]. Figure 6 .

[0136] 4) Set SH(original data X, ν, 10, 5), where v is (-1, -1.5, -2, -2.5, -3). See the results below. Figure 7 .

[0137] Example 4

[0138] 1) A near-infrared spectrometer (PE, integration time 60s, resolution 2cm) was used.-1 Near-infrared spectra of sucrose, fructose, glucose (AR), and raw brown sugar were collected (data interval 1 nm, wavelength range 1000-2500 nm) to obtain the original spectra T, as shown in the figure. Figure 8 As shown.

[0139] 2) Select the fructose spectral data G1, set SH(original data G1,ν,50,3), where ν is order 1-4, with an interval of 1 order. The results are shown in [link to results]. Figure 9 .

[0140] 3) Select the glucose spectral data P1, set SH(original data P1,ν,50,3), where ν is order 1.1-1.5 with an interval of 0.1. The results are shown in [the original data set]. Figure 10 .exist Figure 10 middle Figure 10 (a) shows the results of normalizing the fractional derivatives for all data. Figure 10 (b) is a magnified view of a portion of the image.

[0141] 4) Randomly set SH(original data T, -1.8, 60, 5), and see the integral processing results for the four sugars. Figure 11 .

[0142] The results of processing simulated Gaussian peaks, sin(x) functions, and actual measured spectral data show that the SH method can perform differentiation and integration of arbitrary-order calculus, including integer and fractional orders, meeting the requirements of robustness, accuracy, and convenience.

[0143] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions conceived without inventive effort should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope defined in the claims.

Claims

1. An arbitrary-order calculus method for spectra, comprising the following steps: Step 1: Spectral acquisition to obtain the raw spectrum to be processed. And collect sample spectra; Step 2: Using the SH operator, calculate the arbitrary-order calculus of the spectrum. 1) The measured spectral signal contains both signal and noise, i.e., the measured value = signal + noise; that is, in the frequency domain, the envelope of the measured value is the sum of the signal envelope and the noise envelope. M = S + N (1) For white noise, N is a constant n, so S is the product of the known M and the response function K; S = M·K (2) So Since the real number field is discontinuous, it can be simplified by simulating a hyperbolic function; Redefining K as SH=asech(bω) p (4) Formula (4) includes three parameters: amplitude a, width b for stopping frequency, and peak steepness exponent p. This function combines the constant n into p, which also means that even for the same signal, the response function corresponding to different levels of noise is different, and the numerical setting can be omitted. Equation (5) separates the signal in the frequency domain and then performs an inverse Fourier transform to transfer it to the time domain; Using convolution for time-domain measurements k is the noise reduction convolution kernel. So s=m*k (8) According to the properties of convolution, the νth derivative of the signal is: s (v) =m*k (v) (9) This property applies to real numbers ν; that is, any order of calculus can be implemented on the convolution kernel with sufficient accuracy. The SH function mentioned earlier has three parameters, of which the amplitude parameter 'a' is not mandatory; to maintain the same magnitude before and after convolution, the kernel 'k' has a necessary condition, namely: Therefore, only b and p are necessary in the Sech function; invSH(b,p)=real(iFFT(self(bω) p ))(11) k = invSH / sum(invSH)(12) k can undertake various tasks of calculus of any order; A SH operator (S = SHfd(Mraw,ν,b,p,ON)) capable of performing arbitrary derivative calculus calculations was designed. It is an executable script for MATLAB. The sub-function Sechpf generates a time-domain convolution kernel based on the simulated frequency domain response Sech and the sub-function GLfd, which is used for GL fractional-order calculations. 2) Input the original spectrum and set the parameters. The calculation S = SHfd(Mraw, ν, b, p, ON) has 5 input terms, of which Mraw and the derivative order ν are required. a) The first parameter, Mraw, is the data sequence to be processed, and the second parameter is the calculus order ν. When only the original data Mraw and the derivative order ν are input, the default is that there is no signal filtering and denoising, and it is just a calculator with arbitrary derivative order. b) The third input parameter is the peak width b, which is mainly adjusted based on the signal-to-noise ratio. c) The setting of the fourth parameter p determines the filter's adaptability to various frequencies, similar to the polynomial order in SG; this parameter does not need to be an integer, nor does it restrict the width setting; by fine-tuning this parameter, the true value can be more closely approximated, which greatly increases the flexibility of the method. d) Any value can be entered in the fifth parameter, the purpose of which is to activate the drawing interface used for tuning; 3) Run the program with the parameters set, and output the processed result S.

2. The method according to claim 1, characterized in that, In step one, the sample spectrum includes visible, ultraviolet, fluorescence, Raman, near-infrared, and mid-infrared spectra.

3. The method according to claim 1, characterized in that, In step 2), the order is not necessarily a positive integer; the differential is positive, the integral is negative, and the fractional order is a non-integer.

Citation Information

Patent Citations

  • Variable-order fractional calculus time domain analysis method based on polynomial fitting

    CN107300854A

  • Inverse matrix iterative deconvolution method for spectral resolution enhancement

    WO2020228047A1