CN103678551A - Large-scale medical image retrieval method based on random sparse codes - Google Patents

Abstract

Provided is a large-scale medical image retrieval method based on random sparse codes. The method includes the following steps that firstly, images to be retrieved are input; secondly, the random sparse codes of all images are solved; thirdly, an autocorrelation matrix is computed; fourthly, the solution of an optimization problem is carried out, low dimensional embedding of all the images is solved through characteristic spectral factorization; fifthly, in a low dimensional space, the Euclidean distance between the images to be retrieved and database images is computed, then ascending sort is carried out on the images, and the Euclidean distance is in negative correlation with similarity; sixthly, the images subjected to the ascending sort are output. The large-scale medical image retrieval method is high in retrieval accuracy and retrieval speed by the application of relative comparison constraints and the random sparse codes.

Description

A kind of extensive medical image search method based on Random sparseness coding

Technical field

The present invention relates to a kind of extensive medical image search method.

Background technology

Image retrieval technologies refers to the technology of from the image data base retrieval image list similar to input picture.Existing technology comprises three Main Branches: text based image retrieval technologies, CBIR technology and consider text and the technology of content.Subjectivity tendency and semantic limitation when the limitation of text based technology shows text marking.Content-based retrieval technology is the mainstream technology of current research, but has some challenge: (1), for the image in difference segmentation field, does not have blanket method; (2) feature selecting is one and is difficult to for a long time the key issue of resolving; (3) select which kind of method for measuring similarity; (4), for the large data of large nuber of images, whether image indexing system can respond fast.In these problems, the most insoluble is feature selecting problem and similarity measurement problem, is also the topic receiving publicity for a long time.The technology that considers text and content has had more quantity of information, but the solution degree of problem is still existed to the space of raising.

When the defect of prior art is to carry out feature selecting mainly with field experience, as adopt grey level histogram as feature, adopt texture as feature etc., the measure of similarity is also generally by rule of thumb, is difficult to find enough reasons, versatility is poor.

Summary of the invention

In order to overcome lower, the slow deficiency of retrieval precision of existing medical image search method, the present invention, by using " relatively constraint " and Random sparseness coding, provides a kind of extensive medical image search method based on Random sparseness coding that improves retrieval precision, promotes retrieval rate.

The technical solution adopted for the present invention to solve the technical problems is:

An extensive medical image search method based on Random sparseness coding, said method comprising the steps of:

Step 1: input image to be retrieved;

Step 2: the Random sparseness coding that solves all images;

All images are carried out to pre-service, the large gray scale two dimensional image such as obtain;

By the gray-scale value of every width gray scale two dimensional image, regard a column vector as, with T, represent the matrix of the capable n row of d, wherein the gray-scale value of all pixels of piece image is shown in every list, wherein the sum of n representative image.The basis matrix of expressing the capable m row of d with B, a base vector is shown in every list, what each base vector was expressed is the ash value of all pixels of the random image of selecting of a width; The code coefficient matrix of expressing the capable n row of m with C, it is parameter to be solved; Random sparseness express to wish to minimize || and T-BC||, solve C, and to require values most in the every column element in Matrix C is 0 under the prerequisite of known T and B;

Step 3: calculate autocorrelation matrix;

After trying to achieve Random sparseness encoder matrix C, order

return model to turn in sparse coding matrix:

square formation

be called autocorrelation matrix, described autocorrelation matrix is expressed as the similarity relation between base vector;

Step 4: separate optimization problem, the low-dimensional that solves all images by spectroscopic eigenvalue analysis embeds, and detailed process is as follows:

(4.1) the following optimization problem of model:

Wherein: A=S-SS,

q is used for depositing relative restraint information, and u represents that the one dimension of all view data embeds, and α is parameter;

(4.2), based on method of Lagrange multipliers and KKT theorem, the optimization problem of formula (2) is converted into:

$\mathrm{Au}-\mathrm{\λ}\stackrel{\‾}{Q}u-\mathrm{\μSu}=0---\left(3\right)$

$u_T\stackrel{\‾}{Q}u\≥\mathrm{\α}---\left(4\right)$

u ^TSu＝1 （5）

1 ^TSu＝0 （6）

λ≥0 （7）

$\mathrm{\λ}(u^T\stackrel{\‾}{Q}u-\mathrm{\α})=0---\left(8\right)$

(4.3) order-μ/λ=β, formula (3) is converted into and solves following secular equation:

$\mathrm{Au}=\mathrm{\λ}(\stackrel{\‾}{Q}-\mathrm{\βS})u---\left(9\right)$

Formula (9) solution vector u* has m, and the satisfied relatively solution vector of constraint, also must meet the feasible constraint condition of formula (5), (6), (7), condition (6) is the orthogonality constraint of broad sense, condition (7) is non-negativity constraint, and wherein λ represents the eigenwert of equation (9), and eigenwert must be non-negative, condition (5) takes to take advantage of the method for a factor to meet, to some proper vector u _i, be multiplied by

just can meet formula (5);

(4.4) the solution vector set that meets formula (5) (6) (7) is denoted as { u _i, and according to u _i ^tau _ivalue ascending sort, get at most before k solution vectors, obtain { u _i} ^(k);

(4.5) to { u _i} ^(k)in each vector, calculate

{ the v obtaining _i} ^(k)the k dimension that is called original digital image data embeds, and is expressed in matrix as:

k is input parameter;

Step 5: in lower dimensional space, calculate the Euclidean distance of image to be retrieved and database images, according to the image apart from ascending sort database, described Euclidean distance and similarity negative correlation;

Step 6: the image after output ascending sort.

Further, in described step 2, code coefficient calculates by following expression:

$C_{\mathrm{ij}}=\frac{K(t_j,b_i)}{{\mathrm{\Σ}}_{i\∈\mathrm{rNB}\left(j\right)}K(t_j,b_i)}---\left(1\right)$

Wherein, K (t _j, b _i) representing cosine function, i ∈ rNB (j) represents base vector b _iat image t _jin r nearest neighbours;

Further again, in described step (4.1), the sampling process of relative restraint matrix Q: set two square formation P and Hs large with Q etc., the initial value design of its all elements is 0, in expert's mark, given a pair of relatively constraint: the similarity between image i and image j is larger than the similarity between image i and image k, sets P _ij← P _ij+ 1, P _ji← P _ji+ 1, H _ik← H _ik+ 1, H _ki← H _ki+ 1, last Q is provided by P-H.

Further, in described step (4.3), the value of β must be less than secular equation

eigenvalue of maximum, when setting the value of input parameter β, need first to matrix

carry out Eigenvalues Decomposition.

In described step 5, first solve image to be retrieved and the Euclidean distance of all database images in k dimension space, according to the image in the big or small ascending order array data storehouse of distance, the time complexity of calculating distance is O(n), during sequence, can adopt heapsort method, safeguard a heap that comprises k node, total time complexity is O(nlog (k)).

Beneficial effect of the present invention is mainly manifested in: (1) feasibility is strong.For field image data base, only need to mark a small amount of " relative restraint information ", do not need a large amount of expert's labor capacity.(2) cost-benefit ratio is high.A small amount of artificial mark can bring increasing substantially in retrieval precision.(3) retrieval rate is fast.(4) strong to the retractility of data scale.The time and spatial complexity of whole method is lower, and the database can the adjustment based on parameter method being stretched more than 1,000,000 width images, and its retrieval response time be second grade.

Accompanying drawing explanation

Fig. 1 is the process flow diagram of the extensive medical image search method based on Random sparseness coding.

Fig. 2 is the schematic diagram of 189 width nuclear magnetic resonance images, wherein, (a) is 80 width images, is (b) 88 width images, is (c) 21 width images.

Fig. 3 is the schematic diagram of image to be retrieved.

Fig. 4 selects at random 30 width images as base from 189 width images, the schematic diagram of the random image of selecting.

Fig. 5 is that all elements of matrix S is illustrated as visual schematic diagram.

Fig. 6 introduces 20 pairs of relatively schematic diagram of constraint, wherein, shown in Fig. 6 (a), is i width image and j width image, shown in Fig. 6 (b), is h width image and j width image.

Fig. 7 is the schematic diagram of the top 10 output image of output.

Embodiment

Below in conjunction with accompanying drawing, the invention will be further described.

With reference to Fig. 1～Fig. 7, a kind of extensive medical image search method based on Random sparseness coding, said method comprising the steps of:

Step 1: input image to be retrieved

Step 2: the Random sparseness coding that solves all images

A given image data set, the method of Random sparseness coding refers to all images (is denoted as, T) basis matrix that is expressed as a random little image subset formation of selecting is (referred to as base, be denoted as: the method for the linear combination B), and the coefficient vector of linear combination (is denoted as: be C) sparse.Can it should be noted that the method is not pursued utilizes base B and the coefficient vector C reconstruct original image, but wishes to catch the high-order structural information in image set.For preciseness, spy does following explanation.The no matter image in database or image to be retrieved, we think that it is all etc. large, if do not wait greatly, with the window of fixed size, as mask, carry out pre-service.In addition, suppose that all images are all gray level images, if not gray level image first will be done gray processing and process.

By the gray-scale value of every width two dimensional image, regard a column vector as, with T, represent the matrix of the capable n row of d, wherein the gray-scale value of all pixels of piece image is shown in every list, wherein the sum of n representative image.The basis matrix of expressing the capable m row of d with B, a base vector is shown in every list, what each base vector was expressed is the ash value of all pixels of the random image of selecting of a width.The code coefficient matrix of expressing the capable n row of m with C, it waits to solve.Random sparseness express to wish to minimize || and T-BC||, solve C, and to require values most in the every column element in Matrix C is 0 under the prerequisite of known T and B.The method of standard can be paid huge time overhead when solving this problem, and what the present invention adopted is that a kind of heuristic very efficiently solves this problem.Code coefficient calculates by following expression:

$C_{\mathrm{ij}}=\frac{K(t_j,b_i)}{{\mathrm{\Σ}}_{i\∈\mathrm{rNB}\left(j\right)}K(t_j,b_i)}---\left(1\right)$

Wherein, K (t _j, b _i) representing cosine function, i ∈ rNB (j) represents base vector b _iat image t _jin the nearest neighbours of r in, in the present invention, parameter r gets 4.

Step 3: calculate autocorrelation matrix

After trying to achieve Random sparseness encoder matrix C, make D _ii=Σ _jc _ij, return model to turn in sparse coding matrix:

square formation

be called autocorrelation matrix, its meaning may be interpreted as the similarity relation between paired base vector.

Step 4: separate optimization problem, the low-dimensional that solves all images by spectroscopic eigenvalue analysis embeds

(4.1) the following optimization problem of model:

Wherein: A=S-SS,

q is used for depositing relative restraint information, and u represents that the one dimension of all view data embeds.α is parameter, and its value is larger, requires the more satisfied relatively constraint of optimization problem, and subsequent content can be mentioned how value of parameter.

Below introduce the obtaining value method of relative restraint matrix Q.Set two square formation P and Hs large with Q etc., the initial value design of its all elements is 0.In expert's mark, given a pair of relatively constraint: the similarity between image i and image j is larger than the similarity between image i and image k, sets P _ij← P _ij+ 1, P _ji← P _ji+ 1, H _ik← H _ik+ 1, H _ki← H _ki+ 1.Last Q is provided by P-H.

(4.2), based on method of Lagrange multipliers and KKT theorem, the optimization problem of formula (2) can be converted into:

$\mathrm{Au}-\mathrm{\λ}\stackrel{\‾}{Q}u-\mathrm{\μSu}=0---\left(3\right)$

$u_T\stackrel{\‾}{Q}u\≥\mathrm{\α}---\left(4\right)$

u ^TSu＝1 （5）

1 ^TSu＝0 （6）

λ≥0 （7）

$\mathrm{\λ}(u^T\stackrel{\‾}{Q}u-\mathrm{\α})=0---\left(8\right)$

(4.3) order-μ/λ=β, formula (3) can be converted into and solve following secular equation:

$\mathrm{Au}=\mathrm{\λ}(\stackrel{\‾}{Q}-\mathrm{\βS})u---\left(9\right)$

Formula (9) solution vector u* has m, and the satisfied relatively solution vector of constraint, also must meet the feasible constraint conditions such as formula (5), (6), (7), condition (6) is the orthogonality constraint of broad sense, condition (7) is non-negativity constraint, wherein λ represents the eigenwert of equation (9), and eigenwert must be non-negative.Condition (5) can take to take advantage of the method for a factor to meet, to some proper vector u _i, be multiplied by

just can meet formula (5).

It is worth mentioning that, although comprise input parameter α in above-mentioned optimization problem, do not need to arrange specially α, only input parameter β need be set.This is because can derive according to formula (3), (4), (5), (8) lower bound that β is α.In order to make this optimization problem have feasible solution, can prove, the value of β must be less than secular equation

eigenvalue of maximum.Therefore, when setting the value of input parameter β, need first matrix Q to be carried out to Eigenvalues Decomposition.

(4.4) the solution vector set that meets formula (5) (6) (7) is denoted as { u _i, and according to u _i ^tau _ivalue ascending sort, get at most before k solution vectors, obtain { u _i} ^(k).

(4.5) to { u _i} ^(k)in each vector, calculate

{ the v obtaining _i} ^(k)the k dimension that is called original digital image data embeds, and is expressed in matrix as: V ∈ R ^{k * n}.The k is here input parameter.Step 5: in lower dimensional space, calculate the Euclidean distance of image to be retrieved and database images, ascending sort, described Euclidean distance and similarity negative correlation.

First solve image to be retrieved and the Euclidean distance of all database images in k dimension space, the image according in the big or small ascending order array data storehouse of distance, also can be understood as according to similarity ascending sort.The time complexity of calculating distance is O(n), during sequence, can adopt heapsort method, safeguard a heap that comprises k node.Total time complexity is O(nlog (k)).

Step 6: output ascending sort result, according to and image to be retrieved between similarity descending sort output database in image.

Example: Figure 2 shows that 189 width nuclear magnetic resonance images, each picture consists of the gray-scale value of 160*160, supposes these image construction image data bases.Image shown in Fig. 3 is image to be retrieved.

(1) the hind computation time is with the variation tendency of image data base scale

By the method for adding random noise to this 189 width image, 1890 width, 18900 width, 189000 width, 1890000 width images have been produced respectively.Object is wish to show that algorithm is along with how the time that the increase computing of data scale expends changes.By adopting random 400 width images as base, can obtain the backstage time consumption trend of listing as table 1.Result shows that, for the database that comprises 1,890,000 images, the hind computation time is about 10 hours.The hind computation time is mainly expended on the Random sparseness coding of step 2.It is worth emphasizing that: the result of these hind computations should be stored in disk with document form, if database has renewal, only need to carry out incremental computations.

The table 1 hind computation time used

(2) retrieval expends time in the variation tendency of database scale

On the basis of hind computation, when a picture to be retrieved of input, expend its retrieval time aspect 3: the Random sparseness coding that calculates picture to be retrieved; Spectral factorization solves low-dimensional and embeds; Solve the distance between image the Output rusults that sorts in lower dimensional space.Suppose that the front k of only output opens the diagnostic result of the most similar picture qualifying picture, adopt the fastest sort method, only need to set up and safeguard a Binary Heap that comprises k node, the execution time is O (nlog (k)) magnitude.Table 2 is listed time and the relation between Database size, the wherein k=20 that retrieval one pictures need to expend.The time that result shows to retrieve 20 similar pictures of a pictures from the database of 1,890,000 images is 4.26 seconds, can meet the demand of practical application completely.

A pictures required time of table 2 retrieval

Image library scale (width)	1890	18900	189000	1890000
					Retrieval time (second)	0.038	0.073	0.4576	4.26

(3) state for example the details that realizes of each step below, data set shown in Fig. 2 for this part.

Step 1: input one secondary image to be retrieved, i.e. image shown in Fig. 3.

Step 2: the Random sparseness coding that solves all images.

(1) from 189 width images, the random 30 width images of selecting are as base, and the random image of selecting as shown in Figure 4.

(2) calculate all images to the similarity of this 30 width image.

It is 25600 vector that the image of every width 160*160 is reorganized into length, and the similarity of two width different images adopts cosine distance metric, calculates two included angle cosine values between vector.All images are obtained to the similarity of 30 width base images, obtained the encoder matrix C of a 30*190.

(3) encoder matrix rarefaction.Each row for C, retain r=4 maximum similarity value, and its residual value all gets 0.

Step 3: calculate autocorrelation matrix

By each column element of C divided by row and, then calculate D _ii=Σ _jc _ijwith

make S represent the autocorrelation matrix of 30*30,

for the convenience of showing, all elements of matrix S is illustrated as the image shown in Fig. 5.

Step 4: separate secular equation, the low dimension that solves all images embeds.

(1) introduce 20 pairs of relatively constraints.Be an example as shown in Figure 6, presentation graphs 6(a) two images in are more similar than two images in Fig. 6 (b).Supposing shown in Fig. 6 (a) it is i width image and j width image, is h width image and j width image shown in Fig. 6 (b), sets P _ijand P _jiadd 1, H _jhand H _hjadd 1, other retrain similar, finally according to Q=P-H calculate relatively constraint matrix Q and

.

(2) setup parameter β=β ₀* γ _max, γ wherein _maxfor Generalized Characteristic Equation

eigenvalue of maximum.Here β ₀be set as 0.5.Secular equation solves the eig function that can utilize Matlab.

(3) secular equation of solution formula (9), obtains series of features value and proper vector.

(4) each proper vector is multiplied by a constant factor, makes it meet formula (5)

(5) will meet the proper vector u of formula (5), (6), (7) _iaccording to u _i ^tau _ibig or small ascending sort.

Make k=10.In proper vector after sequence, get front 10 matrix U that form 10 row 30 row, further obtain

so 10 dimensions that each row of V are former images embed.

Step 5: in 10 dimension spaces, solve the Euclidean distance between the embedding that is embedded into other each images of former image.

Step 6: arrange image output according to the big or small ascending order of distance.Figure 7 shows that front 20 images of output, according to similarity from left to right, reduce successively from top to bottom.Can find out that this technology has higher precision.

Claims (5)

1. an extensive medical image search method of encoding based on Random sparseness, is characterized in that: said method comprising the steps of:

Step 1: input image to be retrieved;

Step 2: the Random sparseness coding that solves all images

All images are carried out to pre-service, the large gray scale two dimensional image such as obtain;

By the gray-scale value of every width gray scale two dimensional image, regard a column vector as, the matrix that represents the capable n row of d with T, wherein the gray-scale value of all pixels of piece image is shown in every list, the sum of n representative image wherein, the basis matrix of expressing the capable m row of d with B, a base vector is shown in every list, what each base vector was expressed is the ash value of all pixels of the random image of selecting of a width; The code coefficient matrix of expressing the capable n row of m with C, it is parameter to be solved; Random sparseness express to wish to minimize || and T mono-BC||, solve C, and to require values most in the every column element in Matrix C is 0 under the prerequisite of known T and B;

Step 3: calculate autocorrelation matrix

After trying to achieve Random sparseness encoder matrix C, make D _ii=Σ _jc _ij, return model to turn in sparse coding matrix:

square formation be called autocorrelation matrix, described autocorrelation matrix is expressed as the similarity relation between base vector;

Step 4: separate optimization problem, the low-dimensional that solves all images by spectroscopic eigenvalue analysis embeds, and detailed process is as follows:

(4.1) the following optimization problem of model:

Wherein: A=S-SS,

q is used for depositing relative restraint information, and u represents that the one dimension of all view data embeds, and α is parameter;

(4.2), based on method of Lagrange multipliers and KKT theorem, the optimization problem of formula (2) is converted into:

$\mathrm{Au}-\mathrm{\λ}\stackrel{\‾}{Q}u-\mathrm{\μSu}=0---\left(3\right)$

$u_T\stackrel{\‾}{Q}u\≥\mathrm{\α}---\left(4\right)$

u ^TSu＝1 （5）

1 ^TSu＝0 （6）

λ≥0 （7）

$\mathrm{\λ}(u^T\stackrel{\‾}{Q}u-\mathrm{\α})=0---\left(8\right)$

(4.3) order-μ/λ=β, formula (3) is converted into and solves following secular equation:

$\mathrm{Au}=\mathrm{\λ}(\stackrel{\‾}{Q}-\mathrm{\βS})u---\left(9\right)$

Formula (9) solution vector u* has m, and the satisfied relatively solution vector of constraint, also must meet the feasible constraint condition of formula (5), (6), (7), condition (6) is the orthogonality constraint of broad sense, condition (7) is non-negativity constraint, and wherein λ represents the eigenwert of equation (9), and eigenwert must be non-negative, condition (5) takes to take advantage of the method for a factor to meet, to some proper vector u _i, be multiplied by

just can meet formula (5);

(4.4) the solution vector set that meets formula (5) (6) (7) is denoted as { u _i, and according to u _i ^tau _ivalue ascending sort, get at most before k solution vectors, obtain { u _i} ^(k);

(4.5) to { u _i} ^(k)in each vector, calculate

{ the v obtaining _i} ^(k)the k dimension that is called original digital image data embeds, and is expressed in matrix as: V ∈ R ^{k * n}, k is input parameter;

Step 5: in lower dimensional space, calculate the Euclidean distance of image to be retrieved and database images, the image in ascending sort database, described Euclidean distance and similarity negative correlation;

Step 6: the result of output ascending sort.

2. a kind of extensive medical image search method based on Random sparseness coding as claimed in claim 1, is characterized in that: in described step 2, code coefficient calculates by following expression:

$C_{\mathrm{ij}}=\frac{K(t_j,b_i)}{{\mathrm{\Σ}}_{i\∈\mathrm{rNB}\left(j\right)}K(t_j,b_i)}---\left(1\right)$

Wherein, K (t _j, b _i) representing cosine function, i ∈ rNB (j) represents base vector b _iat image t _jin r nearest neighbours.

3. a kind of extensive medical image search method based on Random sparseness coding as claimed in claim 1 or 2, it is characterized in that: in described step (4.1), the sampling process of relative restraint matrix Q: set two square formation P and Hs large with Q etc., the initial value design of its all elements is 0, in expert's mark, given a pair of relatively constraint: the similarity between image i and image j is larger than the similarity between image i and image k, sets P _ij← P _ij+ 1, P _ji← P _ji+ 1, H _ik← H _ik+ 1, H _ki← H _ki+ 1, last Q is provided by P-H.

4. a kind of extensive medical image search method based on Random sparseness coding as claimed in claim 3, is characterized in that: in described step (4.3), the value of β must be less than secular equation

eigenvalue of maximum, when setting the value of input parameter β, need to first separate this secular equation and obtain eigenvalue of maximum.

5. a kind of extensive medical image search method based on Random sparseness coding as claimed in claim 1 or 2, it is characterized in that: in described step 5, first solve image to be retrieved and the Euclidean distance of all database images in k dimension space, according to the image in the big or small ascending order array data storehouse of distance, the time complexity of calculating distance is O(n), during sequence, can adopt heapsort method, safeguard a heap that comprises k node, total time complexity is O(nlog (k)).

Images