CN102968762A - Polyethylene glycol terephthalate (PET) reconstruction method based on sparsification and Poisson model - Google Patents

Abstract

The invention discloses a polyethylene glycol terephthalate (PET) reconstruction method based on a sparsification and Poisson model. The PET reconstruction method includes: firstly, acquiring projection data, determining image size range and pixel range, and calculating a system probability matrix; obtaining an initial image through a filtered back projection (FBP) traditional algorithm, and using an obtained log-likelihood function as a reconstruction recovery item; using a wavelet transformation and discrete cosine transformation (DCT) mixed base and weighting as sparse regularization constraint, decomposing an objective function by using a split Bregman method to obtain two sub-problems, using a first sub-problem as a sparse regularization problem under a Gaussian model, and using linear Bregman interation for solving; using a second sub-problem as a Poisson denoising problem, and using a close operator method for solving; operating the last one iterative formula, completing an integrated iteration, obtaining a reconstructed image, and serving as an initial value of the next iteration.

Description

A kind of PET method for reconstructing based on rarefaction and Poisson model

Technical field

The present invention relates to positron emission computer tomography field, relate in particular to a kind of PET method for reconstructing based on rarefaction and Poisson model.

Background technology

Positron emission computer tomography (Positron Emission Tomography, PET) be the more advanced clinical examination image technology of the field of nuclear medicine, be the at present only New video technology that can show at live body biomolecule metabolism, acceptor and neurotransmitter activity, now be widely used in diagnosis and the aspects such as antidiastole, state of an illness judgement, therapeutic evaluation, organ function research and new drug development of various diseases.

The process of PET imaging is by injection or takes radiopharmaceutical that the photon number that gathers the positron annihilation generation obtains data for projection.But data acquisition longer duration, the mass storage data of detector rings record is also brought difficulty for follow-up data storage and image reconstruction, often will carry out data by means of computer cluster processes, increased use cost, data acquisition time is longer, injury to human body is also larger, if but the sampling visual angle is inadequate, and pseudo-shadow will appear in the image of reconstruction; Affected by detector efficiency, it is large that the data for projection that detects is disturbed by poisson noise, the not high major defect that becomes the PET system imaging of reconstructed image quality.

In the middle of PET reconstruction field, being considered to comparatively reasonably, the observation data model is Poisson distribution and approximate Gaussian distribution model.The statistical property of portraying observation data with Poisson distribution is the model of closing to reality situation the most, Poisson model has become one of well-known classical model in the image processing field in the reconstruction, and maximum (maximum-likelihood expectation-maximization, MLEM) algorithm of maximum likelihood expectation of following this model to put forward also is one of classic algorithm of PET image reconstruction.Maximum likelihood estimates to be adaptive to the statistical property of poisson noise, obtained using widely, but the increase along with iterations, noise can be exaggerated, these people have been introduced regularization term to be retrained, the sparse canonical of image turns to image reconstruction new regularization thought is provided, what mostly use at present is the function of single base, such as discrete cosine transform, but single base can't carry out the most effective expression to the heterogeneity image, and the image sparse representation method also is the key factor that affects picture quality.All be under the situation of Gauss's least square, to set up objective function mostly in addition, finding the solution of the method fine fitted Gaussian noise of energy of least square, but polluted by quantum noise, quantum noise is obeyed the Poisson distribution Statistical Principles but not additive noise.

Summary of the invention

The object of the present invention is to provide a kind of PET method for reconstructing based on rarefaction and Poisson model, reduce poisson noise to the impact of imaging, effectively improve the quality of PET reconstructed image.

A kind of PET method for reconstructing based on rarefaction and Poisson model comprises following step:

1) obtains data for projection y and the projection probability matrix A of system, y=(y by the PET imaging system ₁, y ₂..., y _M) ^TThe data for projection that expression detects, y ₁, y ₂..., y _MM the data for projection that expression PET detects;

2) to step 1) in data for projection y carry out FBP and rebuild, obtain initial reconstructed image, and definite gradation of image scope and size requirement;

3) utilize described step 1) in data for projection y and the projection probability matrix A of system, set up objective function

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, u ₁, u ₂..., u _NN pixel value of image after expression is rebuild, H (u) are that the image of Poisson likelihood recovers item, and J (u) is sparse regularization term, and λ is regularization parameter; And the reconstructed image behind the iterative and after being optimized.

Described step 3) in:

$\underset{u}{\arg \min }H\left(u\right)=\underset{u}{\arg \min }\underset{i}{\mathrm{\Σ}}({\left[\mathrm{Au}\right]}_i-y_i\ln {\left[\mathrm{Au}\right]}_i)$

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, A={a _IjBe M * N system projection probability matrix, a _IjThe probability that the photon that expression is launched from j pixel is received detector by i, y _iBe y ₁, y ₂..., y _MIn i data, i is the right sequence number of detector, N represents the number of pixels of reconstructed image, M is the number of data for projection.

Described sparse regularization term is:

$J\left(u\right)=\mathrm{\α}{\left|\right|W_{\mathrm{\ψ}}^k\mathrm{\Ψu}\left|\right|}_1+(1-\mathrm{\α}){\left|\right|W_{\mathrm{\φ}}^k\mathrm{\Φu}\left|\right|}_1$

In the formula: ψ represents the discrete cosine transform operator matrix, and Φ represents the wavelet transformation operator matrix,

Be the weighting of described ψ and Φ, u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, α represents the weight of ψ and Φ, k represents iterations.

Use division Bregman method to decompose described step 3) in objective function and obtain the first subproblem, the second subproblem and an iterative relation formula;

Described the first subproblem is: $u^{k+1}=\underset{u}{\arg \min }\mathrm{\λJ}\left(u\right)+\frac{1}{2\mathrm{\γ}}{\left|\right|\mathrm{Au}-b^k+p^k\left|\right|}_2^2$

Described the second subproblem is: $b^{k+1}=\underset{b}{\arg \min }\underset{i}{\mathrm{\Σ}}(b_i-y_i\ln b_i)+\frac{1}{2\mathrm{\γ}}$

Described iterative relation formula is: p ^K+1=p ^k+ (Au ^K+1-bk ^{+ 1})

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, λ is regularization parameter, and J (u) is sparse regularization term, and A is system's projection probability matrix, and b represents auxiliary variable, b _iI the element of expression b, k represents iterations, γ represents relaxation parameter, y _iBe y ₁, y ₂..., y _MIn i data, p represents the Bregman parameter.

Simultaneously, establish initial value u ⁰=fbp (y), b ⁰=Au ⁰, p ⁰=0.

Described first subproblem is converted to sparse Regularization Problem under the Gauss model, and uses linear Bregman iterative, obtain u ^K+1

Described the second subproblem is converted to the Poisson Denoising Problems, and utilizes the method close on operator to find the solution to obtain b ^K+1

Utilize described u ^K+1And b ^K+1Obtain p by the iterative relation formula ^K+1

With described u ^K+1, b ^K+1And p ^K+1As the initial value of next iteration, and carry out loop iteration to finishing to obtain u ^*, and the reconstructed image after the formation optimization.

Compared with prior art, the invention provides a kind of PET method for reconstructing based on rarefaction and Poisson model, use the Bayesian MAP framework, negative logarithm Poisson likelihood function is as recovering item, the repetition weighted sum mixed base of application image strengthens the sparse property of signal, set up the PET image reconstruction objective function of the poisson noise model of sparse regularization, use optimized algorithm division Bregman iteration to find the solution, compare with algorithms of different by experiment, the result shows that this algorithm can improve the quality of reestablishment imaging under low projection angle, reduce poisson noise to the advantages such as impact of imaging.

Description of drawings:

Fig. 1 is the general flow chart that the present invention is based on the PET method for reconstructing of rarefaction and Poisson model.

Fig. 2 is the sub-process figure of first subproblem of the present invention.

Fig. 3 is the Zubal model image of testing formation method.

Fig. 4 is the sinogram of 30 sampling angles.

Fig. 5 is the sinogram of 90 sampling angles.

Fig. 6 does not carry out the reconstruction figure (90 sampling angles) that full variation is adjusted among the present invention.

Fig. 7 is the reconstruction figure (90 sampling angles) that is not weighted rarefaction among the present invention.

Fig. 8 is the figure (90 sampling angles) after the present invention rebuilds.

Fig. 9 is three kinds of algorithm reconstructions of 30 sampling angles figure: (a) the MLEM algorithm is rebuild figure; (b) the PIDAL algorithm is rebuild figure; (c) reconstruction figure of the present invention.

Figure 10 is three kinds of algorithm reconstructions of 90 sampling angles figure: (a) the MLEM algorithm is rebuild figure; (b) the PIDAL algorithm is rebuild figure; (c) reconstruction figure of the present invention.

Figure 11 is the CORR curve map that PIDAL algorithm and the present invention rebuild under 30 sampling angles.

Figure 12 is the SNR curve map that PIDAL algorithm and the present invention rebuild under 30 sampling angles.

Figure 13 is the CORR curve map that PIDAL algorithm and the present invention rebuild under 90 sampling angles.

Figure 14 is the SNR curve map that PIDAL algorithm and the present invention rebuild under 90 sampling angles.

Figure 15 is the reconstruction figure of three kinds of algorithms of 30 sampling angles and the 50th row pixel value figure of model image.

Figure 16 is the reconstruction figure of three kinds of algorithms of 90 sampling angles and the 50th row pixel value figure of model image.

Embodiment:

For improving the quality of PET imaging, reduce noise and keep the edge information partial result, we have invented a kind of PET method for reconstructing based on mixed base and the sparse regularization of weighting, utilize Poisson model to reduce noise, and as shown in Figure 1, the main algorithm step is as follows:

(1) obtains data for projection y by the PET imaging system, computing system projection probability matrix A, y=(y ₁, y ₂..., y _M) ^TThe data for projection that expression detects, y ₁, y ₂..., y _MM the data for projection that expression PET detects.

(2) data in the step (1) are carried out FBP and rebuild, obtain initial reconstructed image, determine gradation of image scope and size requirement.

(3) set up log-likelihood function as the restoration project scalar functions of rebuilding:

$u^{*}=\underset{u}{\arg \min }H\left(u\right)$

$=\underset{u}{\arg \min }\underset{i}{\mathrm{\Σ}}({\left[\mathrm{Au}\right]}_i-y_i\ln {\left[\mathrm{Au}\right]}_i)$

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, A={a _IjBe M * N system projection probability matrix, a _IjThe probability that the photon that expression is launched from j pixel is received detector by i, y _iBe y ₁, y ₂..., y _MIn i data, i is the right sequence number of detector, N represents the number of pixels of reconstructed image, M is the number of data for projection.

Because the data for projection that detects is subjected to the interference of poisson noise, utilizes wavelet transformation, dct transform mixed base and weighting are as sparse regularization constraint item:

$J\left(u\right)=\mathrm{\α}{\left|\right|W_{\mathrm{\ψ}}^k\mathrm{\Ψu}\left|\right|}_1+(1-\mathrm{\α}){\left|\right|W_{\mathrm{\φ}}^k\mathrm{\Φu}\left|\right|}_1$

In the formula: ψ represents the discrete cosine transform operator matrix, and Φ represents the wavelet transformation operator matrix,

The weighting of expression ψ and these two kinds of bases of Φ, α represents the weight of base;

Obtaining final objective function is:

$u^{*}=\underset{u}{\arg \min }H\left(u\right)+\mathrm{\λJ}\left(u\right)$

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, u ₁, u ₂..., u _NN pixel value of image after expression is rebuild, H (u) are that the image of Poisson likelihood recovers item, and J (u) is sparse regularization term, and λ is regularization parameter;

Use division Bregman method that objective function is decomposed two sub-problems and the iterative relation formula of obtaining:

$\left\{\begin{array}{c}{u}^{k+1}=\underset{u}{\arg \min }\mathrm{\λJ}\left(u\right)+\frac{1}{2\mathrm{\γ}}{\left|\right|\mathrm{Au}-b^k+p^k\left|\right|}_2^2\\ b^{k+1}=\underset{b}{\arg \min }\underset{i}{\mathrm{\Σ}}(b_i-y_i\ln b_i)+\frac{1}{2\mathrm{\γ}}{\left|\right|b-{\mathrm{Au}}^{k+1}-p^k\left|\right|}_2^2\\ p^{k+1}=p^k+({\mathrm{Au}}^{k+1}-b^{k+1})\end{array}\right.$

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, λ is regularization parameter, and J (u) is sparse regularization term, and A is system's projection probability matrix, and b represents auxiliary variable, b _iI the element of expression b, k represents iterations, γ represents relaxation parameter, y _iBe y ₁, y ₂..., y _MIn i data, p represents the Bregman parameter, and the p pattern will be used vector matrix.

u ^K+1Be first subproblem, i.e. the regular terms subproblem; b ^K+1Be the second sub-problems, i.e. fidelity item subproblem; p ^K+1Expression iterative relation formula can direct solution.γ＞0 wherein, the complexity of whole algorithm mainly determines by the first subproblem and the second subproblem, and significantly simplifies finding the solution of former problem.

And make initial value u ⁰=fbp (y), b ⁰=Au ⁰, p ⁰=0.

(4) first subproblem is used as sparse Regularization Problem under the Gauss model, uses linear Bregman iterative, with reference to process flow diagram 2, obtain iterative formula; Carry out sub-iteration initialization, make t=b ^k-p ^k, set initial v ₀=0, sub-iterations l, the weighting weight of base

With $W_{\mathrm{\ψ}}^0=1.$

(5) enter sub-iteration, upgrade v _L+1=v _l+ A ^T(t-Au ^k), v is carried out the total variation adjustment:

${v^{\′}}_{l+1}=v_{l+1}-{\mathrm{\β}}_l\frac{\∂\mathrm{TV}\left(v_{l+1}\right)}{{\∂v}_{i,j}}$

v _{I, j}The element of expression v, TV represents total variation, β _lBe step factor, add the total variation adjustment and be in order to utilize the sparse characteristic of image border that it can further improve the image reconstruction quality, wherein the partial derivative of formula can be written as:

$\frac{\∂\mathrm{TV}\left(v\right)}{{\∂v}_{i,j}}=\frac{v_{i,j}-v_{i-1,j}}{\sqrt{{(v_{i,j}-v_{i-1,j})}^2+{(v_{i-1,j+1}-v_{i-1,j})}^2+{\mathrm{\ϵ}}^2}}$

$+\frac{v_{i,j}-v_{i,j-1}}{\sqrt{{(v_{i+1,j-1}-v_{i,j-1})}^2+{(v_{i,j}-v_{i,j-1})}^2+{\mathrm{\ϵ}}^2}}$

$-\frac{v_{i+1,j}+v_{i,j+1}-{2v}_{i,j}}{\sqrt{{(v_{i+,j}-v_{i,j})}^2+{(v_{i,j+1}-v_{i,j})}^2+{\mathrm{\ϵ}}^2}}$

Introduce a very little parameter ε and be in case ask behind the local derviation as infinitely great, the parameter ε of introduction should be less than or equal to 1% v maximal value, the ε value conference smoothly fall the edge.Carry out carrying out soft threshold values processing at wavelet basis and DCT base again after the total variation adjustment:

$\left\{\begin{array}{c}{\mathrm{<figure-callout\; id="1"\; label="u"\; filenames="HDA00002296912800021.png,HDA00002296912800023.png"\; state="\{\{state\}\}">u</figure-callout>}}_{1,l+1}={\mathrm{\&Psi;}}^{-1}\mathrm{shrink}(\mathrm{\&Psi;}{\mathrm{\&delta;v}}_{l+1},{\mathrm{\&delta;\&lambda;w}}_{\mathrm{\&psi;}}^l)\\ {\mathrm{<figure-callout\; id="2"\; label="u"\; filenames="HDA00002296912800062.png"\; state="\{\{state\}\}">u</figure-callout>}}_{2,l+1}={\mathrm{\&Phi;}}^{-1}\mathrm{shrink}(\mathrm{\&Phi;}{\mathrm{\&delta;v}}_{l+1},{\mathrm{\&delta;\&lambda;w}}_{\mathrm{\&phi;}}^l)\end{array}\right.$

In the formula:

β _lBe step factor, ψ ^-1Expression DCT inverse transformation operator, Φ ^-1Expression wavelet inverse transformation operator matrix, shrink are soft threshold values operations, that is:

$\mathrm{shrink}(y,a)=\frac{y}{\left|y\right|}*\max \left\{\right|y|-a,0\}$

And obtain u according to the weight reconstruct of base _L+1:

u _l+1＝αu _1，l+1+(1-α)u _2，l+1。

α is the weight of base.The weighting of two kinds of bases of final updating

With

Can obtain after finishing whole sub-iteration

The product of vector in the algorithm, the wavelet transformation inverse transformation, matrix differentiate difference can realize fast, can be good at obtaining the result.

(6) the second sub-problems is used as the Poisson Denoising Problems, utilizes the method for closing on operator directly to find the solution, close on operator and be defined as:

${\mathrm{prox}}_{\mathrm{\&gamma;f}}:H\&RightArrow;H:x\&RightArrow;\underset{y\&Element;H}{\arg \min }f\left(y\right)+\frac{1}{2\mathrm{\&gamma;}}{\left|\right|y-x\left|\right|}_2^2$

Prox _{γ f}Be the unique solution of above-mentioned conversion, in this problem, can make

Then the second sub-problems can be expressed as:

b ^k+1＝prox _γf(Au ^k+1+p ^k)

According to the definition of contiguous operator, objective function is separable about b, can be converted into One Dimension Optimization Problems and find the solution, and tries to achieve analytical expression and is:

$b_i^{k+1}=\frac{1}{2}\{{[{\mathrm{Au}}^{k+1}+p^k]}_i-\mathrm{\&gamma;}+\sqrt{{{({[{\mathrm{Au}}^{k+1}+p^k]}_i-\mathrm{\&gamma;})}_i}^2+4\mathrm{\&gamma;}{y}_i}\}$

Obtain iterative formula.

(7) last iterative problem of computing: p ^K+1=p ^k+ (Au ^K+1-b ^K+1), finish once complete iteration, the image after obtaining rebuilding, and as the next iteration initial value, until after all iteration is finished, the image u after obtaining rebuilding ^*

The present invention uses the Bayesian MAP framework, negative logarithm Poisson likelihood function is as recovering item, the repetition weighted sum mixed base of application image strengthens the sparse property of signal, set up the PET image reconstruction objective function of the poisson noise model of sparse regularization, use optimized algorithm division Bregman iteration to find the solution, make the precision of imaging obtain certain raising, effectively reduced noise.Fig. 6 uses the present invention but does not carry out the reconstruction figure that full variation is adjusted under 90 sampling angles, Fig. 7 uses the present invention but the reconstruction figure that is not weighted rarefaction under 90 sampling angles, Fig. 8 is the reconstruction figure of the present invention under 90 sampling angles, contrasting this three width of cloth figure can find to adopt full variation to adjust, can carry out smoothly image, edge effect has also obtained maintenance, and the weighting rarefaction can be good at rarefaction, rebuild better effects if, improve the precision of imaging.

Can find out that by Fig. 9, Figure 10 the present invention has better image reconstruction effect, can be under the condition of owing to sample (limited angle) noise of effectively drawing up, compare with the PIDAL algorithm with the MLEM algorithm and to have obvious advantage.Comparison diagram 9 and Figure 10 detect the increase of photon number along with the increase of sampling angle in the simulation process, and the image reconstruction effect has good improvement.But in the PET of reality detection process, the acquisition of sampled data depends on equipment to detection efficiency and the processing power of photon.The reconstruct effect all has some improvement under the sampling condition owing can to find out the present invention.

The simulation results of the present invention

It is the Zubal model of 128 * 128 pixels that emulation experiment of the present invention adopts size, and as shown in Figure 3, the pixel span of template image is between 0 to 5.In the simulation of PET detection data, experimental group is 3 * 10 for the total number of light photons that detects ⁵, be 30 sampling angles, the control experiment group is 9 * 10 for detecting photon number ⁵, be 90 sampling angles.Detection data is obeyed Poisson distribution, and wherein analogue delay is counted r at random _iAccount for 10% of total photon number, set simultaneously correction coefficient c _iThe log-normal stochastic variable that equals 0.3 for obeying standard deviation, the data for projection parameter of simulation is 130 * 30,130 * 90 two kinds, 30,90 projecting directions have namely evenly distributed in 0 to 180 degree scope, each direction has 130 detectors pair, and Fig. 4 and Fig. 5 represent respectively the sinogram under these two kinds of samplings.The linear Bregman iterations that we choose inside in the experiment is 25 times, and outside division Bregman iterations is 40 times, and experiment parameter is by experience manual adjustments, sparse regularization parameter λ=0.02, relaxation parameter γ=0.2, weight parameter α=0.5.The Matlab image reconstruction kit that provides by Fessler etc., all programs are all at Pentium (R) Dual-Core 2.5GHz processor, use Matlab 7.6 to realize on the PC of 2GB internal memory, with existing classical MLEM algorithm, compare with the PIDAL algorithm, each algorithm all is based on this image reconstruction kit and realizes.

Verify by experiment the reconstruction algorithm validity that proposes, relatively rebuild the effect of recovering, judge the quality of each class methods, we adopt following two kinds of criterions: with related coefficient (correlation coefficient, CORR) estimate as the spatial similarity between reconstructed image and the original image and signal to noise ratio (S/N ratio) (signal to noise ratio, SNR) as the standard of judging.Definition CORR is:

$\mathrm{CORR}=\frac{\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}(x_j-u_x)({\mathrm{xtrue}}_j-u_{\mathrm{xtrue}})}{\sqrt{\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{(x_j-u_x)}^2\underset{j=1}{\overset{J}{\mathrm{\&Sigma;}}}{({\mathrm{xtrue}}_j-u_{\mathrm{xtrue}})}^2}}$

X in the formula _jImage pixel value after expression is rebuild, xtrue _jBe the pixel value of realistic model, u _xBe the average of the image pixel after rebuilding, u _XtrueAverage for the image pixel after rebuilding.The larger expression reconstruction of related coefficient effect is better.

Definition SNR is:

$\mathrm{SNR}=10\lg \left\{\frac{\underset{j=1}{\overset{m}{\mathrm{\&Sigma;}}}{(x_j-u_{\mathrm{xtrue}})}^2}{\underset{j=1}{\overset{m}{\mathrm{\&Sigma;}}}{(x_j-{\mathrm{xtrue}}_j)}^2}\right\}$

The quality of the larger explanation image reconstruction of SNR is higher.

Figure 11 and Figure 12 have shown CORR and the SNR curve map of PIDAL algorithm and imaging of the present invention under 30 sampling angles; Figure 13 and Figure 14 have shown CORR and the SNR curve map of PIDAL algorithm and imaging of the present invention under 90 sampling angles.As can be seen from the figure along with the increase image reconstruction effect of the number of times of iteration is improved gradually, but be tending towards at last convergence, algorithm of the present invention all is higher than the PIDAL algorithm on two quality assessment parameters.Figure 15 has shown the 50th road wheel profile of the image after algorithms of different under 30 sampling angles is rebuild and the degree of closeness of original template outline line, and Figure 16 has shown the 50th road wheel profile of the image after algorithms of different under 90 sampling angles is rebuild and the degree of closeness of original template outline line.On Figure 15 and Figure 16, can find out, compare algorithm of the present invention with other algorithm and can play level and smooth effect at regional area better near former figure, rebuild better effects if.

Claims (9)

1. the PET method for reconstructing based on rarefaction and Poisson model is characterized in that, comprises following step:

1) obtain data for projection y and the projection probability matrix A of system by the PET imaging system, and data for projection y=(y ₁, y ₂..., y _M) ^T, y ₁, y ₂..., y _MM the data for projection that expression PET detects;

2) to step 1) in data for projection y carry out FBP and rebuild, obtain initial reconstructed image, and definite gradation of image scope and size requirement;

3) utilize described step 1) in data for projection y and the projection probability matrix A of system, set up objective function In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, u ₁, u ₂..., u _NN pixel value of image after expression is rebuild, H (u) are that the image of Poisson likelihood recovers item, and J (u) is sparse regularization term, and λ is regularization parameter; And the reconstructed image behind the iterative and after being optimized.

2. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 1 is characterized in that described step 3) in:

$\underset{u}{\arg \min }H\left(u\right)=\underset{u}{\arg \min }\underset{i}{\mathrm{\&Sigma;}}({\left[\mathrm{Au}\right]}_i-y_i\ln {\left[\mathrm{Au}\right]}_i)$

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, A={a _IjBe M * N system projection probability matrix, a _IjThe probability that the photon that expression is launched from j pixel is received detector by i, y _iBe y ₁, y ₂..., y _MIn i data, i also is the right sequence number of detector, N represents the number of pixels of reconstructed image, M is the number of data for projection.

3. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 2 is characterized in that, described sparse regularization term is:

$J\left(u\right)=\mathrm{\&alpha;}{\left|\right|W_{\mathrm{\&psi;}}^k\mathrm{\&Psi;u}\left|\right|}_1+(1-\mathrm{\&alpha;}){\left|\right|W_{\mathrm{\&phi;}}^k\mathrm{\&Phi;u}\left|\right|}_1$

In the formula: ψ represents the discrete cosine transform operator matrix, and Φ represents the wavelet transformation operator matrix,

Be the weighting of described ψ and Φ, u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, α represents the weight of ψ and Φ, k represents iterations.

4. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 3 is characterized in that, uses division Bregman method to decompose described step 3) in objective function and obtain the first subproblem, the second subproblem and an iterative relation formula;

Described the first subproblem is: $u^{k+1}=\underset{u}{\arg \min }\mathrm{\&lambda;J}\left(u\right)+\frac{1}{2\mathrm{\&gamma;}}{\left|\right|\mathrm{Au}-b^k+p^k\left|\right|}_2^2$

Described the second subproblem is: $b^{k+1}=\underset{b}{\arg \min }\underset{i}{\mathrm{\&Sigma;}}(b_i-y_i\ln b_i)+\frac{1}{2\mathrm{\&gamma;}}$

Described iterative relation formula is: p ^K+1=p ^k+ (Au ^K+1-b ^K+1)

In the formula: u=(u ₁, u ₂..., u _N) ^TExpression reconstructed image vector, λ is regularization parameter, and J (u) is sparse regularization term, and A is system's projection probability matrix, and b represents auxiliary variable, b _iI the element of expression b, k represents iterations, γ represents relaxation parameter, y _iBe y ₁, y ₂..., y _MIn i data, p represents the Bregman parameter.

5. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 4 is characterized in that, establishes initial value u ⁰=fbp (y), b ⁰=Au ⁰, p ⁰=0.

6. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 5 is characterized in that, described first subproblem is converted to sparse Regularization Problem under the Gauss model, and uses linear Bregman iterative, obtains u ^K+1

7. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 6 is characterized in that, described the second subproblem is converted to the Poisson Denoising Problems, and utilizes the method close on operator to find the solution to obtain b ^K+1

8. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 7 is characterized in that, utilizes described u ^K+1And b ^K+1Obtain p by the iterative relation formula ^K+1

9. the PET method for reconstructing based on rarefaction and Poisson model as claimed in claim 8 is characterized in that, with described u ^K+1, b ^K+1And p ^K+1As the initial value of next iteration, and carry out loop iteration to finishing to obtain u ^*, and the reconstructed image after the formation optimization.

Images