CN109829942B - An automatic quantification method of retinal vessel diameter in fundus images - Google Patents

Abstract

The invention belongs to the technical field of medical image processing, and discloses a method for automatically quantifying retinal vessel diameters in fundus images; using a two-step positioning method that combines thick and thin positioning to automatically detect optic discs; preprocessing the fundus images, including denoising the fundus images and enhanced processing; use multi-scale analysis method to separate retinal blood vessels from fundus image background, use mathematical morphology to open and reconstruct to eliminate false blood vessels in segmentation; use boundary method to measure blood vessel diameter, and display the results. Based on the characteristics of the fundus image, the present invention first detects the optic disc, then enhances the image contrast, segments the retinal blood vessels, and finally measures the size of the retinal blood vessel diameters, displays the results on the system interface, and provides early prediction and assistance to clinicians. Diagnosis provides reference, saves labor and improves efficiency.

Description

An automatic quantification method of retinal vessel diameter in fundus images

technical field

The invention belongs to the technical field of medical image processing, and in particular relates to a method for automatically quantifying retinal vessel diameters in fundus images.

Background technique

At present, the commonly used existing technology in the industry is as follows: the fundus image of retinal blood vessels is an economical and effective clinical auxiliary method for ophthalmic diseases, and its development has attracted widespread attention. Changes in retinal arteriovenous intersections, arterial branches, venous paths, and morphology are important symptoms of retinal arteriosclerosis. When doctors observe fundus images, they mainly use radiography, projection, ophthalmoscopy, and calibration observation to measure the diameter of retinal vessels. Most of these methods are judged by manual observation, manual measurement and experience, and cannot technically accurately detect and extract retinal blood vessels and quantify the diameter of blood vessels, which will not be conducive to early prediction and auxiliary diagnosis of cardiovascular and cerebrovascular diseases, resulting in the I don't have a clear understanding of my condition. Therefore, manual analysis of fundus images has many shortcomings such as strong subjectivity, low diagnostic efficiency, poor repeatability, long time consumption, high labor intensity, and difficulty in ensuring accuracy and consistency. The invention uses digital image processing technology to measure the diameter of retinal blood vessels, overcomes the shortcomings of large subjective errors and poor repeatability, and has the advantages of objective measurement data, repeatable measurement, and labor cost savings.

Contents of the invention

Aiming at the problems existing in the prior art, the present invention provides a method for automatically quantifying the diameter of retinal blood vessels in fundus images.

The present invention is achieved in this way, a method for automatically quantifying the diameter of retinal vessels in fundus images, the method for automatically quantifying the diameters of retinal vessels in fundus images includes the following steps:

In step one, the optic disc is detected using a two-step positioning method combining coarse and fine positioning.

Optic disc positioning is the basis of retinal blood vessel tracking, macular and lesion extraction in fundus images, and various parameters of optic disc are of great significance for early prediction and clinical auxiliary diagnosis. The shape of the edge of the optic disc is between an ellipse and a circle. The rough positioning regards it as a circle, and the fine positioning detects the real shape.

(1) Coarse positioning of optic disc

According to the gray level distribution of the fundus image, binarization is carried out, and some interference is removed by using mathematical morphology knowledge to obtain the optic disc interval, roughly locate the centroid of the optic disc, and roughly position the optic disc as a circle. Assuming that the centroid of the optic disc is O, find the upper and lower points A, B and the left and right points C and D of the optic disc edge along the horizontal and vertical directions of point O, respectively, a total of four points, use the least square method to fit the circle to find the centroid, and then obtain the sub-image of the optic disc.

The coordinates (x _o , y _o ) of the centroid O are calculated as follows,

Among them, x _i and y _i are the abscissa and ordinate of the optic disc edge points respectively, and N is the total number of the optic disc edge points.

Under the principle of the smallest sum of squared errors, the least squares estimate of the circle parameter is:

X＝(U ^T U) ^-1 U ^T V＝U ^-1 V (2)

in

(x ₁ , y ₁ ), (x ₂ , y ₂ ), (x ₃ , y ₃ ) are the coordinates of any three points in A, B, C, and D, and the mean values of the circle center coordinates in the four cases are respectively taken as Optic disc centroid coordinates.

(2) Optic disc fine positioning

Using the optic disc brightness information and local blood vessel information to roughly locate the optic disc, the position of the optic disc can be roughly determined, and the sub-image can be intercepted. Using the dynamic contour tracking method based on non-reinitialized variational level sets can easily obtain the best optic disc boundary curve, which has good fault tolerance and robustness, and realizes precise positioning of the optic disc.

The image energy function can be expressed as

为符号距离函数。

其中d是点(x,y)到曲线的最短距离，

表惩罚项，它确保水平集函数为符号距离函数，μ为惩罚项权重，常系数λ和v，λ>0，δ(■)是Dirac函数，停止函数

H(■)是Heaviside函数。In the formula, Ω represents the image domain,

is a signed distance function.

where d is the shortest distance from the point (x,y) to the curve,

Table penalty item, which ensures that the level set function is a signed distance function, μ is the weight of the penalty item, constant coefficients λ and v, λ>0, δ(■) is the Dirac function, and the stop function

H(■) is a Heaviside function.

Step 2: Preprocessing the fundus image, performing denoising and enhancement processing, using a multi-scale method to separate the retinal blood vessels from the background of the fundus image, and eliminating the pseudo blood vessels that exist in the segmentation.

Distortion occurs in individual images, and the distortion correction is realized by using linear interpolation method and normalized image processing method. In the imaging process of fundus images, due to the influence of non-ideal acquisition conditions and differences in imaging sensors and other characteristics, the collected images often have problems such as noise, uneven illumination, and low image contrast. Therefore, fundus images are preprocessed before retinal vessel segmentation. Preprocessing mainly includes image denoising and contrast enhancement, mathematical morphology processing, etc.

(1) Fundus image denoising and contrast enhancement based on Contourlet transform. In the Contourlet domain, the soft threshold denoising method is used to suppress the noise; the non-linear enhancement operator is used to enhance the transformed sub-band coefficients in each bandpass direction, and the contrast-enhanced image is obtained through inverse transformation.

The nonlinear enhancement operator is as follows:

E(xi _,j )＝xi _,j .sign(xi _,j ).tanh(bu).(1+c.exp(-uu) (4)

in,

(2) Retinal vessel segmentation based on multi-scale analysis linear tracking

The process of segmenting retinal blood vessels by using multi-scale analysis method is the process of extracting blood vessels of different regions and widths at different scales. First select a part of the seed pixels from the image and start tracking. During the tracking process, a higher confidence is given to the pixels that meet the conditions. After the tracking is completed, the confidence matrix of all pixels is quantified to obtain the initial blood vessel.

① Selection of initial seed for linear tracking

Tracing the vessel path of the initial seed is defined as:

V _s ＝{(x,y:T _low <I(x,y)<T _high ) (5)

V _s is the selected seed sequence, I(x, y) represents the gray level of the image pixel at (x, y), T _low and T _high represent the minimum and highest gray value of the image, through the size of the blood vessel area in the image It is estimated that the points with pixels lower than T _low belong to the darker area of retinal image blood vessels, and the points with pixels higher than T _high belong to the brighter area of retinal image blood vessels, and the confidence of these pixel points is higher.

During the linear tracking process, the confidence degree of the pixel belonging to the blood vessel is stored in the sequence _Cw . If the confidence degree of a certain pixel in the confidence matrix is high, the probability that it belongs to the blood vessel is relatively high. Initially, the elements of all scales in the confidence matrix are set to 0.

C _w (x, y): = 0 (6)

②Initialization of linear tracking

The initialization process of linear tracking of retinal blood vessels can be expressed by formula (7):

k:=1, V _c (k):=V _s (t), C _c :={} (7)

Among them, V _c is the pixel sequence obtained by tracking under the current loop variable t, and C _c is the new tracking pixel sequence.

③ Estimation of new tracking pixels

The currently tracked pixel is the last object entering _Vc , and _Cc is a sequence composed of candidate pixels, and the candidate pixels are the nearest 8 pixels around the current tracking pixel, excluding points belonging to _Vc , such as Formula (8) shows:

C _c = N ₈ (V _c (k)) - V _c (8)

Add the current tracking pixel to _Vc , as shown in formula (9):

k:=k+1, V _c (k):=(x,y) (9)

Then the confidence matrix C _w can be updated as:

C _w (x,y):=C _w (x,y)+1,(x,y):=(x _c ,y _c ) (10)

When all parameter values in formula (9) are less than T, start to search for the next seed pixel point (t:=t+1).

④Multi-scale linear tracking

Repeatedly track all seed points according to a certain scale. The number of scales is determined by the width of the blood vessels to be detected in the retinal image. If the blood vessels in the retinal image have greater variability, the number of corresponding scales will also be larger. After several simulation experiments, the effect is better when the scale W=3,5,7,9,11.

⑤ Initial extraction of retinal blood vessels

Initial estimates of blood vessels were performed using map quantification. The initial blood vessels are constructed from pixels larger than T _C in the confidence matrix, and usually the T _C value is equal to the value of the scale W.

Step 3: measure the diameter of the blood vessel using the boundary method, and display the result on the system interface.

(1) Selection of blood vessel measurement image area

From the segmented retinal blood vessel image, select the region of interest to intercept the blood vessel sub-image for edge detection, and then measure the width of the blood vessel in the sub-image, which can reduce the processing time.

(2) Calculation of vessel diameter width

The direction of the blood vessels in the sub-image should be parallel to the horizontal direction as much as possible. If it is not horizontal, the method of rotating the angle can be used to make it parallel to the horizontal direction as much as possible. The minimum distance method is used to realize the automatic measurement of the diameter and width of retinal vessels.

For the coordinates of the upper and lower edge points in the vertical direction, use the distance formula between two points to calculate, and calculate the distance value as the blood vessel width value, which can be expressed as:

Wherein, w _i is the width value of column i, and (x _i , y _i ) and (x _k , y _k ) are the coordinates of the upper and lower edges of the blood vessel respectively. The distance method has the advantages of wide adaptability, good stability and high precision.

(3) Ruler positioning and pixel spacing calculation

In order to convert the vascular diameter marked by pixels into the actual size of the vascular caliber, the pixel pitch needs to be measured, therefore, the scale must be positioned in the image.

The physical length of the pixel is calculated by the number of pixels between the minimum tick value and the adjacent tick mark. Set the scale in the image of the same size as the intercepted sub-image, for example, the scale is 10mm×10mm, after performing edge detection processing on the scale image, a binary image is obtained, the pixel value of the scale is set to 1, and the pixel value of the background is 0. For the binary scale image, traverse the binary image from top to bottom, and calculate the number of pixels.

决定，p的单位为mm/pixel。Let the size of the fixed ruler be 10mm, and the number of pixels between the scales be Num, then the physical length p of the pixels can be given by

Determined, the unit of p is mm/pixel.

(4) Convert the pixel size of the vessel diameter to the physical size

By multiplying the pixels measured by each marker point with p, the physical size of the diameter of the retinal blood vessels can be known.

Another object of the present invention is to provide an automatic quantification system for realizing the retinal blood vessel diameter measurement method of the fundus image, the system comprising:

The fundus image optic disc detection module is used to realize the detection of the fundus image optic disc;

The preprocessing and blood vessel segmentation module is used to achieve fundus image denoising and contrast enhancement, and to achieve retinal blood vessel segmentation;

The blood vessel width measurement and application module is used to realize retinal blood vessel diameter measurement, early prediction and auxiliary diagnosis for clinicians.

To sum up, the advantages and positive effects of the present invention are: on the basis of detecting the optic disk of the fundus image, the image is preprocessed to improve the contrast of the fundus image. Multi-scale analysis method is used to segment retinal blood vessels, and the boundary method is used to measure the diameter of the vessels, and the results are displayed on the system interface to provide references for early prediction and auxiliary diagnosis of clinicians, saving labor and improving efficiency.

Description of drawings

Fig. 1 is a flowchart of a method for automatically quantifying retinal vessel diameters in fundus images provided by an embodiment of the present invention.

2 is a schematic structural diagram of an automatic quantification system for retinal blood vessel diameters in fundus images provided by an embodiment of the present invention;

In the figure: 1. Image preprocessing module; 2. Detection and segmentation module; 3. Quantization and application module.

Fig. 3 is a flowchart for realizing the method for automatically quantifying retinal vessel diameters in fundus images provided by an embodiment of the present invention.

Fig. 4 is a schematic diagram of the effect of the rough positioning process provided by the embodiment of the present invention.

In the figure: a, optic disc after interference removal; b, schematic diagram of centroid positioning; c, subimage of optic disc; d, rough positioning result.

Fig. 5 is a schematic diagram of the fine positioning effect provided by the embodiment of the present invention.

In the figure: a, dynamic tracking process diagram; b, tracking result diagram; c, tracking curve bold; d, fine positioning result.

Fig. 6 is a schematic diagram of the effect of preprocessing and retinal vessel segmentation provided by the embodiment of the present invention.

In the figure: a, the original image; b, the contrast-enhanced image; c, the initial segmentation of blood vessels; d, the blood vessels after removing interference.

Fig. 7 is a schematic diagram of the measurement effect of the retinal blood vessel width measurement process provided by the embodiment of the present invention.

In the figure: a, intercepting blood vessels; b, edge detection; c, removing isolated noise; d, rotating blood vessels; e, establishing a ruler; f, measurement results.

Detailed ways

In order to make the object, technical solution and advantages of the present invention more clear, the present invention will be further described in detail below in conjunction with the examples. It should be understood that the specific embodiments described here are only used to explain the present invention, not to limit the present invention.

Aiming at the problems of manual observation and manual measurement of blood vessel diameter in the prior art, labor intensity is high and time-consuming is long. The invention automatically measures the diameter of the retinal blood vessels through the image processing technology, provides early prediction and reference for auxiliary diagnosis for doctors, reduces the labor intensity of doctors, and improves diagnosis efficiency.

The application principle of the present invention will be described in detail below in conjunction with the accompanying drawings.

As shown in Figure 1, the method for automatically quantifying the diameter of retinal blood vessels in fundus images provided by the embodiment of the present invention includes the following steps:

S101: Automatically detect the optic disc by using a two-step positioning method combining coarse and fine positioning;

S102: Preprocessing the fundus image, separating retinal blood vessels from the background of the fundus image, and removing pseudo blood vessels in the segmentation;

S103: Measure the diameter of the blood vessel by the boundary method, and display the result on the system interface.

In the step S101, the automatic detection of the optic disc location in the optic disc is the basis for retinal blood vessel tracking, macular and lesion extraction in fundus images, and various parameters of the optic disc are of great significance to clinical diagnosis. The shape of the edge of the optic disc is between an ellipse and a circle. The rough positioning regards it as a circle, and the fine positioning detects the real shape.

(1) Coarse positioning of optic disc

According to the gray level distribution of the fundus image, the image is binarized, and the knowledge of mathematical morphology is used to remove part of the interference, obtain the optic disc interval, roughly locate the centroid of the optic disc, and roughly position the optic disc as a circle. Assuming that the centroid of the optic disc is O, find the upper and lower points A, B and the left and right points C and D of the optic disc edge along the horizontal and vertical directions of point O, respectively, a total of four points, use the least square method to fit the circle to find the centroid, and then obtain the sub-image of the optic disc.

The coordinates (x _o , y _o ) of the centroid O are calculated as follows,

Among them, x _i and y _i are the abscissa and ordinate of the optic disc edge points respectively, and N is the total number of the optic disc edge points.

Under the principle of the smallest sum of squared errors, the least squares estimate of the circle parameter is:

X＝(U ^T U) ^-1 U ^T V＝U ^-1 V (2)

in

(x ₁ , y ₁ ), (x ₂ , y ₂ ), (x ₃ , y ₃ ) are the coordinates of any three points in A, B, C, and D, and the mean values of the circle center coordinates in the four cases are respectively taken as Optic disc centroid coordinates.

(2) Optic disc fine positioning

Using the optic disc brightness information and local blood vessel information to roughly locate the optic disc, the position of the optic disc can be roughly determined, and the sub-image can be intercepted. Using the dynamic contour tracking method based on non-reinitialized variational level sets can easily obtain the best optic disc boundary curve, which has good fault tolerance and robustness, and realizes precise positioning of the optic disc.

The image energy function can be expressed as

为符号距离函数。

其中d是点(x,y)到曲线的最短距离，

表惩罚项，它确保水平集函数为符号距离函数，μ为惩罚项权重，常系数λ和v，λ>0，δ(■)是Dirac函数，停止函数

H(■)是Heaviside函数。In the formula, Ω represents the image domain,

is a signed distance function.

where d is the shortest distance from the point (x,y) to the curve,

Table penalty item, which ensures that the level set function is a signed distance function, μ is the weight of the penalty item, constant coefficients λ and v, λ>0, δ(■) is the Dirac function, and the stop function

H(■) is a Heaviside function.

In the step S102, if some images are distorted, linear interpolation method and normalized image processing method are used to realize distortion correction. In the imaging process of fundus images, due to the influence of non-ideal acquisition conditions and differences in imaging sensors and other characteristics, the collected images often have problems such as noise, uneven illumination, and low image contrast. Therefore, fundus images are preprocessed before retinal vessel segmentation. Preprocessing mainly includes image denoising and contrast enhancement, mathematical morphology processing, etc.

(1) Fundus image denoising and contrast enhancement based on Contourlet transform. In the Contourlet domain, the soft threshold denoising method is used to suppress the noise; the non-linear enhancement operator is used to enhance the transformed sub-band coefficients in each bandpass direction, and the contrast-enhanced image is obtained through inverse transformation.

The nonlinear enhancement operator is as follows:

E(xi _,j )＝xi _,j .sign(xi _,j ).tanh(bu).(1+c.exp(-uu) (4)

in,

(2) Retinal vessel segmentation based on multi-scale analysis linear tracking

The process of segmenting retinal blood vessels by using multi-scale analysis method is the process of extracting blood vessels of different regions and widths at different scales. First select a part of the seed pixels from the image and start tracking. During the tracking process, a higher confidence is given to the pixels that meet the conditions. After the tracking is completed, the confidence matrix of all pixels is quantified to obtain the initial blood vessel.

① Selection of initial seed for linear tracking

Tracing the vessel path of the initial seed is defined as:

V _s ＝{(x,y:T _low <I(x,y)<T _high ) (5)

Vs is the selected seed sequence, I(x, y) represents the gray level of the image pixel at (x, y), T _low and T _high are the minimum gray value and the highest gray value of the image, through the blood vessel area in the image Estimated by the size of , the points with pixels lower than T _low belong to the darker area of retinal image blood vessels, and the points with pixels higher than T _high belong to the brighter area of retinal image blood vessels, and the confidence of these pixel points is higher.

During the linear tracking process, the confidence degree of the pixel belonging to the blood vessel is stored in the sequence _Cw . If the confidence degree of a certain pixel in the confidence degree matrix is high, the probability that it belongs to the blood vessel is relatively high. Initially, the elements of all scales in the confidence matrix are set to 0.

C _w (x, y): = 0 (6)

②Initialization of linear tracking

The initialization process of retinal blood vessel linear tracking can be expressed by formula (7):

k:=1, V _c (k):=V _s (t), C _c :={} (7)

Among them, V _c is the pixel sequence obtained by tracking under the current loop variable t, and C _c is the new tracking pixel sequence.

③ Estimation of new tracking pixels

The currently tracked pixel is the last object entering _Vc , and _Cc is a sequence composed of candidate pixels, and the candidate pixels are the nearest 8 pixels around the current tracking pixel, excluding points belonging to _Vc , such as Formula (8) shows:

C _c = N ₈ (V _c (k)) - V _c (8)

Add the current tracking pixel to _Vc , as shown in formula (9):

k:=k+1, V _c (k):=(x,y) (9)

Then the confidence matrix C _w can be updated as:

C _w (x,y):=C _w (x,y)+1,(x,y):=(x _c ,y _c ) (10)

When all parameter values in formula (9) are less than T, start to search for the next seed pixel point (t:=t+1).

④Multi-scale linear tracking

Repeatedly track all seed points according to a certain scale. The number of scales is determined by the width of the blood vessels to be detected in the retinal image. If the blood vessels in the retinal image have greater variability, the number of corresponding scales will also be more. After several simulation experiments, the effect is better when the scale W=3,5,7,9,11.

⑤ Initial extraction of retinal blood vessels

Initial estimates of blood vessels were performed using map quantification. The initial blood vessels are constructed from pixels larger than T _C in the confidence matrix, and usually the value of T _C is equal to the scale W value.

In the step S103, the specific process of measuring the vascular diameter by the boundary method is as follows:

(1) Selection of blood vessel measurement image area

From the segmented retinal blood vessel image, select the region of interest to intercept the blood vessel sub-image for edge detection, and then measure the width of the blood vessel in the sub-image, which can reduce the processing time.

(2) Calculation of vessel diameter width

The direction of blood vessels in the sub-image should be parallel to the horizontal direction as much as possible. If it is not horizontal, the method of rotation angle can be used to make it parallel to the horizontal direction as much as possible. The minimum distance method is used to realize the automatic measurement of the diameter and width of retinal vessels.

For the coordinates of the upper and lower edge points in the vertical direction, use the distance formula between two points to calculate, and calculate the distance value as the blood vessel width value, which can be expressed as:

Wherein, w _i is the width value of column i, (x _i , y _i ) and (x _j , y _j ) are the coordinates of the upper and lower edges of the blood vessel respectively. The distance method has the advantages of wide adaptability, good stability and high precision.

(3) Ruler positioning and pixel spacing calculation

In order to convert the vascular caliber marked by pixels into the actual size of the vascular caliber, the pixel pitch needs to be measured, therefore, the scale must be positioned in the image.

The physical length of the pixel is calculated by the number of pixels between the minimum tick value and the adjacent tick mark. Set the scale in the image of the same size as the intercepted sub-image, the scale is 10mm×10mm, after the edge detection process is performed on the scale image, a binary image is obtained, the pixel value of the scale is set to 1, and the pixel value of the background is 0 . For the binary scale image, traverse the binary image from top to bottom, and calculate the number of pixels.

决定，p的单位为mm/pixel。Let the size of the fixed ruler be 10mm, and the number of pixels between the scales be Num, then the physical length p of the pixels can be given by

Determined, the unit of p is mm/pixel.

(4) Convert the pixel size of the vessel diameter to the physical size

By multiplying the pixels measured by each marker point with p, the physical size of the diameter of the retinal blood vessels can be known.

As shown in FIG. 2 , the system for automatic quantification of retinal vessel diameters in fundus images provided by the embodiment of the present invention includes: an optic disc detection module 1 , a preprocessing and vessel segmentation module 2 , and a vessel width measurement and application module 3 .

The optic disc detection module 1 is used to realize retinal image disc detection.

Preprocessing and vessel segmentation module 2 is used to improve image quality and extract retinal vessels.

The blood vessel width measurement and application module 3 is used to measure the diameter of retinal blood vessels and provide references for clinicians.

The application principle of the present invention will be described in detail below in conjunction with specific embodiments.

Retinal Optic Disc Detection:

As shown in FIG. 4 , the rough positioning process provided by the embodiment of the present invention: find the centroid of the optic disc, and after finding the centroid, roughly position it as a circle, and intercept the sub-image.

As shown in FIG. 5 , the fine positioning process provided by the embodiment of the present invention is: dynamic tracking process diagram, tracking result diagram, tracking curve thickening, and fine positioning result.

Preprocessing and retinal vessel segmentation:

As shown in FIG. 6 , the preprocessing and retinal vessel segmentation process provided by the embodiment of the present invention includes: original image, contrast-enhanced image, initial vessel segmentation, and vessel after interference removal.

Retinal Vessel Width Measurement:

As shown in FIG. 7 , the measurement process of the retinal blood vessel width provided by the embodiment of the present invention includes: intercepting the blood vessel, detecting the edge, removing isolated noise, rotating the blood vessel, establishing a scale, and measuring the result.

The above descriptions are only preferred embodiments of the present invention, and are not intended to limit the present invention. Any modifications, equivalent replacements and improvements made within the spirit and principles of the present invention should be included in the protection of the present invention. within range.

Claims (2)

1. a method for automatic quantification of retinal blood vessel caliber in fundus image, characterized in that, the automatic quantification method for retinal blood vessel caliber in described fundus image comprises the following steps: Step 1, the optic disc is detected by a two-step positioning method combining coarse and fine positioning; Step 2, preprocessing the fundus image, performing denoising and enhancement processing, using a multi-scale method to separate the retinal blood vessels from the background of the fundus image, and eliminating the pseudo blood vessels that exist in the segmentation; Step 3, using the boundary method to measure the diameter of the blood vessel, and display the result on the system interface; In the first step, the positioning of the optic disc in the automatic detection of the optic disc is the basis of retinal blood vessel tracking, macula and lesion extraction in fundus images: (1) Coarse positioning of optic disc Binarize the image according to the gray distribution of the fundus image, use mathematical morphology knowledge to remove part of the interference, obtain the optic disc interval, roughly locate the centroid of the optic disc, and roughly position the optic disc as a circle; assume that the centroid of the optic disc is O, respectively along the O In the horizontal and vertical directions of the points, find the upper and lower points A, B and the left and right points C, D on the edge of the optic disc, a total of four points, use the least square method to fit the circle to find the center of mass, and then obtain the sub-image of the optic disc; The coordinates (x _o , y _o ) of the centroid O are calculated as follows:

Wherein x _i and y _i are the abscissa and ordinate of the optic disc edge points respectively, and N is the total number of the optic disc edge points; Under the principle of the minimum sum of squared errors, the least squares estimation of the circle parameters is: X＝(U ^T U) ^-1 U ^T V＝U ^-1 V (2) in

(x ₁ , y ₁ ), (x ₂ , y ₂ ), (x ₃ , y ₃ ) are the coordinates of any three points in A, B, C, and D, and the mean values of the circle center coordinates in the four cases are respectively taken as Optic disc centroid coordinates; (2) Optic disc fine positioning Use the optic disc brightness information and local blood vessel information to roughly locate the optic disc, and use the dynamic contour tracking method based on the non-reinitialized variational level set to conveniently obtain the best optic disc boundary curve, which has good fault tolerance and robustness, and realizes the optic disc The precise location of: The image energy function is expressed as:

In the formula, Ω represents the image domain,

is a signed distance function;

where d is the shortest distance from the point (x,y) to the curve,

Table penalty item, which ensures that the level set function is a signed distance function, μ is the weight of the penalty item, constant coefficients λ and v, λ>0, δ(■) is the Dirac function, and the stop function

H( ) is the Heaviside function; In the second step, distortion occurs in individual images, and the distortion correction is realized by using the linear interpolation method and the normalized image processing method; during the imaging process of the fundus image, due to the influence of non-ideal acquisition conditions and differences in imaging sensor characteristics, the collected images There are problems of noise, uneven illumination, and low image contrast; preprocessing mainly includes image denoising, contrast enhancement, and mathematical morphology processing; (1) Fundus image denoising and contrast enhancement based on Contourlet transform; in the Contourlet domain, the soft threshold denoising method is used to suppress noise; the non-linear enhancement operator is used to enhance the transformed sub-band coefficients in each band-pass direction, and after inversion transform to obtain a contrast-enhanced image; The nonlinear enhancement operator is as follows: E(xi _,j )＝xi _,j .sign(xi _,j ).tanh(bu).(1+c.exp(-uu) (4) in,

(2) Retinal vessel segmentation based on multi-scale analysis linear tracking The image is expressed in multiple scales, and the multi-scale analysis method is used to segment the retinal blood vessels into blood vessels of different regions and widths at different scales; first, a part of the seed pixels are selected from the image to start tracking, and during the tracking process, the satisfying conditions are selected The pixels are given higher confidence, and after the tracking is completed, the confidence matrix of all pixels is quantified to obtain the initial blood vessel; ① Selection of initial seed for linear tracking Tracing the vessel path for the initial seed is defined as: V _s ＝{(x,y:T _low <I(x,y)<T _high ) (5) V _s is the selected seed sequence, I(x, y) represents the gray level of the image pixel at point (x, y), T _low and T _high are the minimum gray value and the highest gray value of the image, through the image The size of the blood vessel area is estimated. The points with pixels lower than T _low belong to the darker area of retinal blood vessels, and the points with pixels higher than T _high belong to the brighter area of retinal image blood vessels. The confidence of these pixels is higher; In the process of linear tracking, the confidence degree of the pixel belonging to the blood vessel is stored in the sequence _Cw . If the confidence degree of a certain pixel in the confidence matrix is high, the probability that it belongs to the blood vessel is relatively high; initially, the confidence degree The elements of all scales in the matrix are set to 0; C _w (x, y): = 0 (6) ②Initialization of linear tracking The initialization process of linear tracking of retinal blood vessels is expressed by formula (7): k:=1, V _c (k):=V _s (t), C _c :={} (7) Among them, V _c is the pixel sequence obtained by tracking under the current loop variable t, and C _c is the new tracking pixel sequence; ③ Estimation of new tracking pixels The currently tracked pixel is the last object that enters _Vc , and _Cc is a sequence composed of candidate pixels, and the candidate pixels are the nearest 8 pixels around the current tracking pixel, excluding _Vc , as shown in the formula As shown in (8): C _c = N ₈ (V _c (k)) - V _c (8) Add the current tracking pixel to _Vc , as shown in formula (9): k:=k+1, V _c (k):=(x,y) (9) Then the confidence matrix C _w is updated as: C _w (x,y):=C _w (x,y)+1,(x,y):=(x _c ,y _c ) (10) When all parameter values are less than T in the formula (9), start to look for the next seed pixel point (t:=t+1); ④Multi-scale linear tracking Repeatedly track all seed points according to a certain scale. The number of scales is determined by the width of the blood vessels to be detected in the retinal image. If the blood vessels in the retinal image have greater variability, the number of corresponding scales is also relatively large; after multiple simulations For experiments, scales W=3, 5, 7, 9, 11 are more appropriate; ⑤ Initial extraction of retinal blood vessels Use the mapping quantification method to initially estimate the blood vessels; the initial blood vessels are constructed from pixels larger than T _C in the confidence matrix, and the T _C value is equal to the value of the scale W; In the third step, the specific process of measuring the vascular diameter by the boundary method is as follows: (1) Selection of blood vessel measurement image area From the segmented retinal blood vessel image, select the region of interest to intercept the blood vessel sub-image for edge detection, and then measure the width of the blood vessel in the sub-image to reduce the processing time; (2) Calculation of vessel diameter width The direction of the blood vessels in the sub-image should be parallel to the horizontal direction as much as possible. If it is not horizontal, use the method of rotation angle to make it parallel to the horizontal direction as much as possible; use the minimum distance method to realize the automatic measurement of the diameter and width of retinal blood vessels; For the coordinates of the upper and lower edge points in the vertical direction, use the distance formula between two points to calculate, and calculate the distance value as the blood vessel width value, expressed as:

Among them, w _i is the width value of column i, (x _i , y _i ) and (x _j , y _j ) are the coordinates of the upper and lower edges of blood vessels respectively; the distance method has the advantages of wide adaptability, good stability and high precision; (3) Ruler positioning and pixel spacing calculation In order to convert the vascular caliber calibrated by pixels into the actual size of the vascular caliber, it is necessary to measure the pixel spacing, therefore, it is necessary to solve the positioning of the ruler in the image; Calculate the physical length of the pixel by the number of pixels between the minimum scale value and the adjacent scale line; set the scale in the image of the same size as the intercepted sub-image, the scale is 10mm×10mm, and perform edge detection processing on the scale image , what is obtained is a binary image, the pixel value of the scale is set to 1, and the pixel value of the background is 0; the binary image is traversed from top to bottom for the binary scale image, and the number of pixels is calculated; Let the size of the fixed ruler be 10mm, and the number of pixels between the scales be Num, then the physical length p of the pixels is given by

Decision, the unit of p is mm/pixel; (4) Convert the pixel size of the vessel diameter to the physical size The physical size of the diameter of the retinal blood vessels can be known by multiplying the pixel measured by each marker point with p. 2. A fundus image retinal vessel diameter automatic quantification system that realizes the method for automatically quantifying the retinal vessel diameter of the fundus image of claim 1, wherein the automatic quantification system for the retinal vessel diameter of the fundus image comprises: The optic disc detection module is used to realize the optic disc detection of the fundus image; Image preprocessing and segmentation module, used to achieve fundus image denoising and contrast enhancement, retinal vessel segmentation and edge detection; The quantification and application module is used to realize retinal vessel diameter measurement, early clinical prediction and auxiliary diagnosis.

Images