CN105203498A - Near infrared spectrum variable selection method based on LASSO - Google Patents

Abstract

The invention provides a near infrared spectrum variable selection method based on LASSO. The method comprises the following concrete processes of collecting a near infrared spectrum of a sample, and using a conventional method for measuring a concentration vector of a tested ingredient; dividing a data set into a training set and a prediction set by adopting a certain grouping mode; determining the constraint value t of an LASSO method through crossed verification; using a minimum angle regression algorithm for calculating the regression coefficient beta; remaining the position of a wavelength point of the beta being not zero; building a partial least squares regression model between a training set spectrum and the concentration vector by utilizing the training set spectrum corresponding to the remained wavelength, and predicting the concentration of a tested ingredient of a prediction set sample. The method has the advantages that an effective wavelength can be extracted; a quantitative analysis model is simplified; the prediction precision of the model is improved. Compared with an existing variable selection method, the method has the advantages that the speed is high; the repeating performance can be realized; the higher prediction precision can be reached by using fewer variables. The near infrared spectrum variable selection method is applicable to the variable selection of complicated sample near infrared spectrums.

Description

A kind of near infrared spectrum Variable Selection based on LASSO

Technical field

This method invention belongs to the nondestructive analysis technical field in analytical chemistry field, is specifically related to a kind of near infrared spectrum Variable Selection based on LASSO.

Background technology

Near-infrared spectral analysis technology is the technology of high speed development in analytical chemistry field, and it has, and analysis efficiency is high, detection speed fast, without the need to advantages such as sample pretreatments, be widely used in the industry such as food, oil.Modling model between the content or classification of near infrared spectrum and measured matter, can realize the direct qualitative and quantitative analysis of complex material.In near infrared spectrum modeling, a very important problem is exactly there is redundancy wavelength in spectrum.General near infrared spectrum (NIR) comprises hundreds of wavelength variable point, and the character of some of them wavelength and research is incoherent, these uncorrelated wavelength points, can affect model quality, causes its predictive ability to decline.Therefore variables choice is the important content of spectrum modeling analysis always.

Variable Selection conventional in spectral data analysis mainly comprises the method for method based on intelligent optimization algorithm and Corpus--based Method.The former mainly contains simulated annealing (simulatedannealing, SA, see SwierengaH, deGrootPJ, deWeijerAP, DerksenMWJ, BuydensLMC, ImprovementofPLSmodeltransferabilitybyrobustwavelengthse lection, ChemomIntellLabSyst, 1998, 41, 237-248), genetic algorithm (geneticalgorithm, GA, see LeardiR, GonzalezAL, GeneticalgorithmsappliedtofeatureselectioninPLSregressio n:howandwhentousethem, ChemomIntellLabSyst, 1998, 41, 195-207), tabu search (Tabusearch, TS, see HagemanJA, StreppelM, WehrensR, WavelengthselectionwithTabuSearch, JChemometrics, 2003, 17, 427-437), ant group algorithm (antcolonyoptimization, ACO, see ShamsipurM, Zare-ShahabadiV, HemmateenejadB, AkhondM, Antcolonyoptimization:apowerfultoolforwavelengthselectio n, JChemometrics, 2006, 20, 146-157), particle cluster algorithm (particleswarmoptimization, PSO, see XuL, JiangJH, WuHL, ShenGL, YuRQ, Variable-weightedPLS, ChemomIntellLabSyst, 2007, 85, 140-143) etc., these optimized methods exist needs a large amount of parameters, search time is longer and be easily absorbed in the defects such as local optimum.The latter mainly contains without information variable removing method (UninformativeVariableElimination, UVE, see CentnerV, MassartDL, deNoordOE, JongS, VandeginsteBM, SternaC, Eliminationofuninformativevariablesformultivariatecalibr ation.AnalChem, 1996, 68, 3851-3858), Monte Carlo combines without information variable removing method (MonteCarloUninformativeVariableElimination, MCUVE, see CaiWS, LiYK, ShaoXG, Avariableselectionmethodbasedonuninformativevariableelim inationformultivariatecalibrationofnear-infraredspectra, ChemomIntellLabSyst, 2008, 90, 188-194), based on the Variable Selection method (RandomizationTest of randomized test, RT, see XuH, LiuZC, CaiWS, ShaoXG, Awavelengthselectionmethodbasedonrandomizationtestfornea r-infraredspectralanalysis.ChemomIntellLabSyst, 2009, 97, 189-193) etc.UVE method have employed leave one cross validation to obtain variable stability value, the repetitious computing of this process need, but also need introduce with original spectrum comprise the equal random noise variable of variables number, so when data set number is larger, the method counting yield is low, consuming time longer.MCUVE algorithm and RT method all introduce repeatedly modeling technique, the multiple models produced often more than single model can effectively extract from the different aspect of data and different aspects and express the complex relationship between independent variable and dependent variable, are conducive to more rationally, reliably choice variable.But due to the Stochastic choice of each modeling sample, make the operation result of these two kinds of methods there is certain instability, and also more time-consuming when data volume is larger.Therefore, be necessary to further develop novel Variable Selection fast, improve stability and the precision of prediction of model.

Summary of the invention

The object of the invention is for above-mentioned existing problems, a kind of quick, stable Variable Selection is provided.The method, under the absolute value sum of a regression coefficient is less than the condition of a constant, makes residual sum of squares (RSS) minimize, thus more strictly makes some regression coefficient vanishing, and corresponding variable is deleted, realizes variables choice.

Concrete steps are as follows:

(1) m sample to be tested is collected.Setting spectrum parameter, the near infrared spectrum of collecting sample, obtains the spectrum matrix X of sample.Measure the tested component concentration of sample by conventional method, obtain concentration vector y.Adopt certain packet mode that data are divided into training set and forecast set, wherein training set sample is used for Modling model Optimal Parameters, and forecast set sample is used for the predictive ability of testing model.

(2) cross validation is adopted to determine the binding occurrence t of LASSO.T controls the degree of compression, and t is less, and the degree of compression is stronger, and due to this restrictive condition, end product can make some component of regression coefficient β become 0, reaches the object of variables choice.

(3) utilize minimum angle regression algorithm to solve the regression coefficient β of LASSO, preserve the wavelength points position that regression coefficient is not 0.

$\begin{array}{cc}\widehat{\mathrm{\β}}=\underset{\mathrm{\β}\∈R^p}{\arg min}\left\{{(y-X\mathrm{\β})}^T(y-X\mathrm{\β})\right\}& s.t.\underset{t=1}{\overset{p}{\Σ}}\left|{\mathrm{\β}}_t\right|\≤t\end{array}$

Minimum angle regression algorithm process is as follows:

1. Renewal model is selected in variables set (activeset), calculates related coefficient absolute value

${\hat{y}}_0=0;{\hat{c}}_{kj}=x_j^T(y-{\hat{y}}_{k-1});{\hat{C}}_k=max\left\{\right|{\hat{c}}_{kj}\left|\right\}$

Upgrade activesetA (k),

$A\left(k\right)=A(k-1)+\left\{\hat{j}\right\};A\left(0\right)=\mathrm{\φ};\hat{j}=\underset{j\∉A\left(k-1\right)}{\arg min}\left\{\left|{\hat{c}}_{kj}\right|\right\}$

2. minimum angular direction (u is determined _k)

Make X _k=(... s _jx _j) _{j ∈ A (k)}

Wherein $s_j=sign\left\{{\hat{c}}_{kj}\right\},{\mathrm{\ω}}_k=A_k{\left(X_k^{T}X\right)}^{-1}{1}_k,A_k={\left(1_k^T{\left(X_k^{T}X\right)}^{-1}{1}_k\right)}^{-0.5}$

1 _kbe important be 1 vector, its length equals | A|.Calculate minimum angular direction: u _k=X _kω _k3. material calculation

When $j\∉A\left(k\right),$ Order $a_{kj}=x_j^T{u}_k$

If | A|=d, then algorithm stops.

Otherwise ${\widehat{\mathrm{\γ}}}_k={\min }_{j\∉A\left(k\right)}^{+}\{{\hat{C}}_k-{\hat{c}}_{kj}/(A_k-a_{kj}),({\hat{C}}_k+{\hat{c}}_{kj})/(A_k+a_{kj})\}$

4. predicated response

$\widetilde{\mathrm{\γ}}=\underset{{\mathrm{\γ}}_j>0,j\∈A\left(k\right)}{min}\left\{{\mathrm{\γ}}_j\right\},$ Wherein ${\mathrm{\γ}}_j=-{\widehat{\mathrm{\β}}}_j/\left(s_j{\mathrm{\ω}}_{kj}\right);{\widetilde{\mathrm{\γ}}}_1=\mathrm{\∞}$

If ${\widetilde{\mathrm{\&gamma;}}}_k<{\widehat{\mathrm{\&gamma;}}}_k,$ Then ${\hat{y}}_k={\hat{y}}_{k-1}+{\widetilde{\mathrm{\&gamma;}}}_k{u}_k$

As .j ∈ A, ${\widehat{\mathrm{\&beta;}}}_j\&LeftArrow;{\widehat{\mathrm{\&beta;}}}_j+\widetilde{\mathrm{\&gamma;}}{\mathrm{\&omega;}}_{kj}{s}_j,$ Otherwise $\widehat{\mathrm{\&beta;}}=0$

$A(k+1)=A\left(k\right)-\left\{\tilde{j}\right\},$ Wherein $\tilde{j}=\underset{j}{\arg min}\left\{{\mathrm{\&gamma;}}_j\right\}$

${\hat{c}}_{k+1,j}=x_j^T\left(y-{\hat{y}}_k\right),$ And ${\hat{C}}_{k+1}=\underset{j}{max}\left\{\left|{\hat{c}}_{k+1,j}\right|\right\},$ Return and perform step 1..

Otherwise ${\hat{y}}_k={\hat{y}}_{k-1}+{\widehat{\mathrm{\&gamma;}}}_k{u}_k$

As j ∈ A, ${\widehat{\mathrm{\&beta;}}}_j\&LeftArrow;{\widehat{\mathrm{\&beta;}}}_j+{\widehat{\mathrm{\&gamma;}}}_k{\mathrm{\&omega;}}_{kj}{s}_j,$ Otherwise ${\widehat{\mathrm{\&beta;}}}_j=0$ Return and perform step 1..

(4) according to the wavelength points position retained, only retain the corresponding wavelength row of training set spectrum matrix, obtain new spectrum matrix, and set up partial least squares regression (PLS) model with training set sample measured component concentration vector.Wherein PLS model because of subnumber by Monte Carlo Cross-Validation in conjunction with F inspection determine.Utilize this model, measure the levels of the tested composition of forecast set sample.

Compared with existing Variable Selection, the present invention has that travelling speed is fast, choice variable has repeatable advantage, and can reach with less variable number and better predict the outcome.

Accompanying drawing explanation

Fig. 1: the near infrared light spectrogram of tobacco sample

Fig. 2: tobacco near infrared spectrum data training set carries out residual sum of squares (RSS) (SSR) mean value of 1000 cross validations and the variance variation diagram along with normalized binding occurrence t, and wherein vertical line represents t value corresponding to optimization model

Fig. 3: the regression coefficient β that after tobacco near infrared spectrum data training set carries out LASSO variables choice, all variablees are corresponding

Fig. 4: UVE, MCUVE, RT, LASSO tetra-kinds of Variable Selections retain the distribution plan of variable

Fig. 5: the near infrared light spectrogram of sesame oil and soybean oil, Rice oil triple blend sample

Fig. 6: the spectroscopic data training set of sesame oil and soybean oil, Rice oil triple blend sample carries out residual sum of squares (RSS) (SSR) mean value of 1000 cross validations and the variance variation diagram along with normalized binding occurrence t, and wherein vertical line represents t value corresponding to optimization model

Fig. 7: the regression coefficient β that after sesame oil carries out LASSO variables choice with soybean oil, Rice oil triple blend sample light modal data training set, all variablees are corresponding

Fig. 8: UVE, MCUVE, RT, LASSO tetra-kinds of Variable Selections retain the distribution plan of variable

Embodiment

For better understanding the present invention, below in conjunction with embodiment the present invention done and describe in detail further, but the scope of protection of present invention being not limited to the scope that embodiment represents.

Embodiment 1:

The present embodiment is applied to near-infrared spectrum analysis, measures the content of reducing sugar value in tobacco sample.Concrete step is as follows:

(1) gather the near infrared spectrum data of tobacco leaf sample, use BrukerVector22/N near infrared spectrometer (German Bruker optical instrument company) to test 269 tobacco sheet samples in different tobacco leaf producing region.NIR light spectrum wave-number range is 4000 ~ 9000cm ^-1, sampling interval is 4 wave numbers, totally 1296 wavelength points, and the near infrared light spectrogram of sample as shown in Figure 1.In tobacco sample, reducing sugar (ReducingSugar) content adopts AAIII type Continuous Flow Analysis instrument (German BranLuebbe company) to measure according to standard method.Before modeling, tobacco leaf sample is divided into two parts at random, comprises training set and forecast set sample, wherein training set sample is used for Modling model, forecast set sample and is used for the predictive ability of testing model.

(2) cross validation is adopted to determine the binding occurrence t of LASSO.T controls the degree of compression, and t is less, and the degree of compression is stronger, and this restrictive condition makes some component of vectorial β become 0, thus reaches the object of variables choice.The present embodiment training set carry out residual sum of squares (RSS) (SSR) mean value of 1000 cross validations and variance along with normalized binding occurrence t change as shown in Figure 2, wherein vertical line represents t value corresponding to optimization model, is 0.103.

(3) the regression coefficient β of LASSO is solved.Utilize minimum angle regression algorithm to solve the regression coefficient β of LASSO, preserve the wavelength points position that regression coefficient is not 0.

$\begin{array}{cc}\widehat{\mathrm{\&beta;}}=\underset{\mathrm{\&beta;}\&Element;R^p}{\arg min}\left\{{(y-X\mathrm{\&beta;})}^T(y-X\mathrm{\&beta;})\right\}& s.t.\underset{t=1}{\overset{p}{\&Sigma;}}\left|{\mathrm{\&beta;}}_t\right|\&le;t\end{array}$

This embodiment regression coefficient β value that after carrying out LASSO variables choice, all variablees are corresponding as shown in Figure 3.

(4) according to the wavelength points position retained, only retain the corresponding wavelength row of training set spectrum matrix, obtain new spectrum matrix.Spectrum matrix and training set sample measured component concentration vector set up partial least squares regression (PLS) model, wherein PLS model because of subnumber by Monte Carlo Cross-Validation in conjunction with F inspection determine.Utilize this model, measure the levels of the tested composition of forecast set sample.This embodiment determine because of subnumber be 8.

UVE, MCUVE, RT, LASSO tetra-kinds of Variable Selections retain the distribution plan of variable as shown in Figure 4.As can be seen from Figure 4, on the one hand, LASSO is roughly the same with the variable range of other three kinds of method choice, which illustrates the rationality of LASSO method choice variable.On the other hand, the variable number that LASSO selects is more less than other three kinds of Variable Selections, and this embodies the superiority of the method.

In order to compare the effect of four kinds of variables choice further, table 1 gives tobacco near infrared data and does not adopt variables choice and set up the modeling effect of PLS model after adopting variables choice.From data in table, LASSO choice variable only 27 is nearly 1/10th of other three kinds of Variable Selections.Computing time 11.89, although slower than the PLS not carrying out variables choice, obviously faster than other Variable Selection.The RMSEP value that LASSO-PLS modeling obtains is minimum, and R value is maximum, illustrates that the method more can improve the precision of prediction of model.Therefore, LASSO-PLS choice variable number compared with other modeling method is few, and computing time is shorter, and precision of prediction is higher.

The results contrast of the different modeling method of table 1 tobacco near infrared data

Embodiment 2:

The present embodiment is applied to near-infrared spectrum analysis, measures the near infrared spectrum data of sesame oil and soybean oil, Rice oil triple blend.Concrete step is as follows:

(1) the NIR light modal data of sesame oil and soybean oil, Rice oil triple blend sample is gathered, use near infrared spectrometer (TJ270-60, Tianjin Tuopu Equipment Co., Ltd.) carry out near infrared spectrum data measurement, wavelength coverage is 800 ~ 2500nm, sampling interval is 1nm, totally 1701 wavelength points.The near infrared light spectrogram of sample as shown in Figure 5.Sample configures by a certain percentage (soybean oil quality 0.05 ~ 2.5, interval 0.05; Rice oil concentration 0.05 ~ 2.5, interval 0.05).Before modeling, sample is divided into two parts at random, comprises training set and forecast set sample, wherein training set sample is used for Modling model, forecast set sample and is used for the predictive ability of testing model.

(2) cross validation is adopted to determine the binding occurrence t of LASSO.T controls the degree of compression, and t is less, and the degree of compression is stronger, and this restrictive condition makes some component of vectorial β become 0, thus reaches the object of variables choice.This embodiment training set carry out residual sum of squares (RSS) (SSR) mean value of 1000 cross validations and variance along with normalized binding occurrence t variation diagram as shown in Figure 6, wherein vertical line represents t value corresponding to optimization model is 0.254.

(3) the regression coefficient β of LASSO is solved.Utilize minimum angle regression algorithm to solve the regression coefficient β of LASSO, preserve the wavelength points position that regression coefficient is not 0.

$\begin{array}{cc}\widehat{\mathrm{\&beta;}}=\underset{\mathrm{\&beta;}\&Element;R^p}{\arg min}\left\{{(y-X\mathrm{\&beta;})}^T(y-X\mathrm{\&beta;})\right\}& s.t.\underset{t=1}{\overset{p}{\&Sigma;}}\left|{\mathrm{\&beta;}}_t\right|\&le;t\end{array}$

The regression coefficient β value that after this embodiment training set carries out LASSO variables choice, all variablees are corresponding as shown in Figure 7.

(4) according to the wavelength points position retained, only retain the corresponding wavelength row of training set spectrum matrix, obtain new spectrum matrix.Spectrum matrix and training set sample measured component concentration vector set up partial least squares regression (PLS) model, wherein PLS model because of subnumber by Monte Carlo Cross-Validation in conjunction with F inspection determine.Utilize this model, measure the levels of the tested composition of forecast set sample.This embodiment determine because of subnumber be 8.

UVE, MCUVE, RT, LASSO tetra-kinds of Variable Selections retain the distribution plan of variable as shown in Figure 8.As can be seen from Figure 8, LASSO is roughly the same with the variable range of other three kinds of method choice, which illustrates the rationality of LASSO method choice variable.On the other hand, the variable number that LASSO selects is more less than other three kinds of Variable Selections, and this embodies the superiority of the method.

In order to compare the effect of four kinds of variables choice further, table 2 gives sesame oil and soybean oil, Rice oil triple blend near infrared spectrum data do not adopt variables choice and set up the modeling effect of PLS model after adopting variables choice.From data in table, LASSO choice variable only 11, far less than the variable that other Variable Selections are selected.2.48 seconds computing times, obviously faster than other Variable Selection.The RMSEP value that LASSO-PLS modeling obtains is minimum, and R value is maximum.Therefore, LASSO-PLS choice variable number compared with other modeling method is few, and computing time is shorter, and precision of prediction is higher.

The results contrast of the different modeling method of table 2 vegetable oil NIR data

Claims (4)

1., based on a near infrared spectrum Variable Selection of LASSO, it is characterized in that comprising following steps:

1) gather the near infrared spectrum data of measured object sample, measure the measured component concentration content of sample in training set by conventional method, adopt certain packet mode that data are divided into training set and forecast set;

2) the binding occurrence t. of LASSO is determined;

3) minimum angle regression algorithm is utilized to solve the regression coefficient β of LASSO;

4) by training set spectrum Matrix Regression factor beta be not 0 wavelength row set up partial least squares regression (PLS) model with concentration vector, utilize this model, the content of prediction unknown sample composition.

2. a kind of near infrared spectrum Variable Selection based on LASSO according to claim 1, is characterized in that: the described detailed process utilizing minimum angle regression algorithm to solve the regression coefficient β of LASSO is:

1. Renewal model is selected in variables set (activeset), calculates related coefficient absolute value

${\hat{y}}_0=0;{\hat{c}}_{kj}=x_j^I(y-{\hat{y}}_{k-1});{\hat{C}}_k=\max \left\{\right|{\hat{c}}_{kj}\left|\right\}$

Upgrade activesetA (k)

$A\left(k\right)=A(k-1)+\left\{\hat{j}\right\};$ A(0)＝φ $\hat{j}=\underset{j\&NotElement;A(k-1)}{\arg \min }\left\{\right|{\hat{c}}_{kj}\left|\right\}$

2. minimum angular direction (u is determined _k)

Make X _k=(... s _jx _j...) _{j ∈ A (k)}

Wherein $s_j=\mathrm{si}gn\left\{{\hat{c}}_{kj}\right\},$ ${\mathrm{\&omega;}}_k=A_k{\left(X_k^{I}X\right)}^{-1}{1}_k,$ $A_k={\left(1_k^T{\left(X_k^{T}X\right)}^{-1}{1}_k\right)}^{-0.5}$

1 _kbe important be 1 vector, its length equals | A|

Calculate minimum angular direction: u _k=X _kω _k

3. material calculation

When $j\&NotElement;A\left(k\right),$ Order $a_{kj}=x_j^T{u}_k$

If | A|=d, then algorithm stops

Otherwise ${\widehat{\mathrm{\&gamma;}}}_k={\min }_{j\&NotElement;A\left(k\right)}^{+}\{{\hat{C}}_k-{\hat{c}}_{kj}/(A_k-a_{kj}),({\hat{C}}_k+{\hat{c}}_{kj})/(A_k+a_{kj})\}$

4. predicated response

$\widetilde{\mathrm{\&gamma;}}=\underset{{\mathrm{\&gamma;}}_j>0,j\&Element;A\left(k\right)}{min}\left\{{\mathrm{\&gamma;}}_j\right\},$ Wherein ${\mathrm{\&gamma;}}_j=-{\widehat{\mathrm{\&beta;}}}_j/\left(s_j\mathrm{\&omega;}{kj}_{}\right);$ ${\widetilde{\mathrm{\&gamma;}}}_l=\mathrm{\&infin;}$

If ${\widetilde{\mathrm{\&gamma;}}}_k<{\widehat{\mathrm{\&gamma;}}}_k,$ Then ${\hat{y}}_k={\hat{y}}_{k-1}+{\widehat{\mathrm{\&gamma;}}}_k{u}_k$

As j ∈ A, ${\widehat{\mathrm{\&beta;}}}_j\&LeftArrow;{\widehat{\mathrm{\&beta;}}}_j+\widetilde{\mathrm{\&gamma;}}{\mathrm{\&omega;}}_{\mathrm{kj}}{s}_j,$ Otherwise $\widehat{\mathrm{\&beta;}}=0$

$A(k+1)=A\left(k\right)-\left\{\tilde{j}\right\},$ Wherein $\tilde{j}=\underset{j}{\mathrm{argmin}}\left\{{\mathrm{\&gamma;}}_j\right\}$

${\hat{c}}_{k+1j}=x_j^I(y-\hat{y}{k}_{}),$ And ${\hat{C}}_{k+1}=\underset{j}{\max }\left\{\right|{\hat{c}}_{k+1j}\left|\right\},$ Return and perform step 1.

Otherwise ${\hat{y}}_k={\hat{y}}_{k-1}+{\widehat{\mathrm{\&gamma;}}}_k{u}_k$

As j ∈ A, otherwise return and perform step 1..

3. a kind of near infrared spectrum Variable Selection based on LASSO according to claim 1, it is characterized in that: the defining method of the binding occurrence t. of described LASSO is cross validation, t controls the degree of compression, t is less, the degree of compression is stronger, this restrictive condition makes some component of vectorial β become 0, thus reaches the object of variables choice.

4. a kind of near infrared spectrum Variable Selection based on LASSO according to claim 1, is characterized in that: described PLS model because of subnumber defining method be that Monte Carlo Cross-Validation is checked in conjunction with F.