U.S. patent application number 15/398116 was filed with the patent office on 2017-05-04 for method and system for image analysis.
The applicant listed for this patent is Straxcorp Pty Ltd. Invention is credited to Anne BOHTE, Ali GHASEMZADEH, Eleanor Jean MACKIE, Aloys MBALA, Ego SEEMAN, Roger ZEBAZE.
Application Number | 20170119333 15/398116 |
Document ID | / |
Family ID | 43731867 |
Filed Date | 2017-05-04 |
United States Patent
Application |
20170119333 |
Kind Code |
A1 |
ZEBAZE; Roger ; et
al. |
May 4, 2017 |
METHOD AND SYSTEM FOR IMAGE ANALYSIS
Abstract
A method for quantifying the effects of a treatment of a
subject, comprising: analysing a sample of the subject, the sample
comprising a material of interest and at least one adjacent
material that is generally adjacent to the material of interest,
the material of interest having a generally different density from
the adjacent material. The analysing comprises: (i) automatically
defining a plurality of regions of interest within an image of the
sample that includes a junction between the material of interest
and the adjacent material, each of said regions of interest having
a width of one or more pixels or voxels; (ii) determining
respective profiles of density, intensity or attenuation within the
regions of interest, each of the respective profiles crossing the
junctions between the material of interest and the adjacent
material; (iii) automatically identifying the material of interest
in the profiles; (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and (v) locating the junction in each of the
S-like portions, including: (a) identifying the greatest difference
in values of the respective profile between an adjacent peak and
trough in a segment of the respective profile that includes at
least some of the material of interest and at least some of the
adjacent material; and (b) locating a point of inflexion in said
segment. This method may also be employed in a method of treatment
of a disease in a subject.
Inventors: |
ZEBAZE; Roger; (Coburg,
AU) ; MACKIE; Eleanor Jean; (Fitzroy, AU) ;
BOHTE; Anne; (East Bentleigh, AU) ; MBALA; Aloys;
(Epping, AU) ; GHASEMZADEH; Ali; (Greensborough,
AU) ; SEEMAN; Ego; (Toorak, AU) |
|
Applicant: |
Name |
City |
State |
Country |
Type |
Straxcorp Pty Ltd |
Melbourne |
|
AU |
|
|
Family ID: |
43731867 |
Appl. No.: |
15/398116 |
Filed: |
January 4, 2017 |
Related U.S. Patent Documents
|
|
|
|
|
|
Application
Number |
Filing Date |
Patent Number |
|
|
15191079 |
Jun 23, 2016 |
9566038 |
|
|
15398116 |
|
|
|
|
14720004 |
May 22, 2015 |
9378566 |
|
|
15191079 |
|
|
|
|
13395379 |
May 24, 2012 |
9064320 |
|
|
PCT/AU2010/001181 |
Sep 10, 2010 |
|
|
|
14720004 |
|
|
|
|
Current U.S.
Class: |
1/1 |
Current CPC
Class: |
A61B 6/505 20130101;
A61B 6/5205 20130101; G06T 11/003 20130101; A61B 6/032 20130101;
G06T 7/0012 20130101; G06T 7/12 20170101; G06T 2207/20104 20130101;
G06T 2207/30008 20130101; G06T 2207/10081 20130101; G06T 2207/30096
20130101; G06T 2200/24 20130101; G06T 2207/10088 20130101; G06T
7/70 20170101 |
International
Class: |
A61B 6/00 20060101
A61B006/00; A61B 6/03 20060101 A61B006/03 |
Foreign Application Data
Date |
Code |
Application Number |
Sep 11, 2009 |
AU |
2009904407 |
Claims
1. A method for quantifying the effects of a treatment of a
subject, comprising: analysing a sample of the subject, the sample
comprising a material of interest and at least one adjacent
material that is generally adjacent to the material of interest,
the material of interest having a generally different density from
the adjacent material, the analysing comprising: (i) automatically
defining a plurality of regions of interest within an image of the
sample that includes a junction between the material of interest
and the adjacent material, each of said regions of interest having
a width of one or more pixels or voxels; (ii) determining
respective profiles of density, intensity or attenuation within the
regions of interest, each of the respective profiles crossing the
junctions between the material of interest and the adjacent
material; (iii) automatically identifying the material of interest
in the profiles; (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and (v) locating the junction in each of the
S-like portions, including: (a) identifying the greatest difference
in values of the respective profile between an adjacent peak and
trough in a segment of the respective profile that includes at
least some of the material of interest and at least some of the
adjacent material; and (b) locating a point of inflexion in said
segment.
2. A method as claimed in claim 1, further comprising applying the
treatment to the subject.
3. A method as claimed in claim 2, including performing the
analysing of the sample before, during and/or after the
treatment.
4. A method as claimed in claim 1, wherein the treatment comprises
treatment with a drug that treats vessel calcification, acute
myocardial infarction, fracture or stroke.
5. A method as claimed in claim 1, wherein the treatment comprises
treatment with a drug that treats bone and/or modifies bone
architecture.
6. A method as claimed in claim 5, wherein the treatment comprises
an anabolic therapy.
7. A method as claimed in 5, comprising determining any differences
in effects of different treatments on bone architecture.
8. A method as claimed in claim 5, comprising assessing an effect
of the treatment on bone using at least one result of the analysing
of the sample of bone.
9. A method as claimed in claim 1, further comprising diagnosing a
disease including the analysing of the sample.
10. A method of assessing a structure of bone during a treatment in
a subject, the method comprising: treating the bone; and analysing
a sample of bone of the subject, the sample comprising a material
of interest that includes the structure and at least one adjacent
material that is generally adjacent to the material of interest,
the material of interest having a generally different density from
the adjacent material, the analysing comprising: (i) automatically
defining a plurality of regions of interest within an image of the
sample that includes a junction between the material of interest
and the adjacent material, each of said regions of interest having
a width of one or more pixels or voxels; (ii) determining
respective profiles of density, intensity or attenuation within the
regions of interest, each of the respective profiles crossing the
junctions between the material of interest and the adjacent
material; (iii) automatically identifying the material of interest
in the profiles; (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and (v) locating the junction in each of the
S-like portions, including: (a) identifying the greatest difference
in values of the respective profile between an adjacent peak and
trough in a segment of the respective profile that includes at
least some of the material of interest and at least some of the
adjacent material; and (b) locating a point of inflexion in said
segment.
11. A method as claimed in claim 10, comprising assessing an effect
of the treatment using at least one result of the analysing of the
sample of bone.
12. A method as claimed in claim 10, comprising assessing bone
fracture risk of the subject.
13. A method as claimed in claim 10, comprising one or more of: i)
detecting normal bone by analysing the structure of the material of
interest; ii) detecting diseased bone by analysing the structure of
the material of interest; and iii) detecting decayed or abnormal
bone by analysing the structure of the material of interest.
14. A method as claimed in claim 10, comprising detecting decayed
or abnormal bone by analysing the structure of the material of
interest, including quantifying cortical architecture, wherein a
cortical compartment comprises compact-cortex, outer transitional
zone, inner transitional zone or any combination thereof.
15. A method as claimed in claim 14, including: i) quantifying
cortical architecture by determining porosity, cortical thickness,
cortical area, cortical density; or ii) quantifying cortical
architecture by determining porosity, cortical thickness, cortical
area, cortical density and distinguishing cortical porosity
attributable to age-related bone decay from normal porosity.
16. A method as claimed in claim 10, comprising detecting decayed
or abnormal bone by analysing the structure of the material of
interest, including quantifying trabecular architecture.
17. A method as claimed in claim 16, including quantifying
trabecular architecture by determining trabecular density,
trabecular number, trabecular thickness, trabecular number,
trabecular bone volume fraction.
18. A method as claimed in claim 10, comprising detecting decayed
or abnormal bone: i) by analysing the structure of the material of
interest, including combining cortical and trabecular architecture;
ii) by analysing the structure of the material of interest, wherein
the decayed or abnormal bone comprises both an abnormal cortex and
an abnormal trabecular compartment; or iii) by analysing the
structure of the material of interest, including determining any
differences in effects of different treatments on bone
architecture.
19. A method of diagnosing or monitoring a disease or the effects
of a disease in a subject, the method comprising: analysing a
sample of the subject, the sample comprising a material of interest
and at least one adjacent material that is generally adjacent to
the material of interest, the material of interest having a
generally different density from the adjacent material, the
analysing comprising: (i) automatically defining a plurality of
regions of interest within an image of the sample that includes a
junction between the material of interest and the adjacent
material, each of said regions of interest having a width of one or
more pixels or voxels; (ii) determining respective profiles of
density, intensity or attenuation within the regions of interest,
each of the respective profiles crossing the junctions between the
material of interest and the adjacent material; (iii) automatically
identifying the material of interest in the profiles; (iv)
analysing said profiles by automatically defining, for each of said
profiles, a series of sections of said respective profile, each of
said sections comprising at least one S-like portion of the
respective profile, the S-like portion containing the junction
between the material of interest and the adjacent material; and (v)
locating the junction in each of the S-like portions, including:
(a) identifying the greatest difference in values of the respective
profile between an adjacent peak and trough in a segment of the
respective profile that includes at least some of the material of
interest and at least some of the adjacent material; and (b)
locating a point of inflexion in said segment.
20. A method as claimed in claim 19, comprising diagnosing vessel
calcification, acute myocardial infarction, fracture or stroke.
21. A method as claimed in claim 19, comprising monitoring a
disease or the effects of a disease, wherein the subject is being
treated or has been treated with a drug that treats vessel
calcification, acute myocardial infarction, fracture or stroke.
22. A method as claimed in claim 19, comprising diagnosing
fracture-vulnerable bone.
23. A method as claimed in claim 19, comprising one or more of: i)
detecting normal bone by analysing the structure of the material of
interest; ii) detecting diseased bone by analysing the structure of
the material of interest; and iii) detecting decayed or abnormal
bone by analysing the structure of the material of interest.
24. A method as claimed in claim 19, comprising detecting decayed
or abnormal bone by analysing the structure of the material of
interest, including quantifying cortical architecture, wherein a
cortical compartment comprises compact-cortex, outer transitional
zone, inner transitional zone or any combination thereof.
25. A method as claimed in claim 24, including: i) quantifying
cortical architecture by determining porosity, cortical thickness,
cortical area and cortical density; or ii) quantifying cortical
architecture by determining porosity, cortical thickness, cortical
area and cortical density and distinguishing cortical porosity
attributable to age-related bone decay from normal porosity.
26. A method as claimed in claim 19, comprising detecting decayed
or abnormal bone: i) by analysing the structure of the material of
interest, including quantifying trabecular architecture; ii) by
analysing the structure of the material of interest, including
quantifying trabecular architecture by determining trabecular
density, trabecular number, trabecular thickness, trabecular number
and trabecular bone volume fraction; iii) by analysing the
structure of the material of interest, including combining cortical
and trabecular architecture; or iv) by analysing the structure of
the material of interest, wherein the abnormal bone comprises an
abnormal cortex and an abnormal trabecular compartment.
27. A method as claimed in claim 19, comprising assessing an effect
of a treatment of the subject using at least one result of the
analysing of the sample.
28. A method of treatment of a disease in a subject, the method
comprising: quantifying the effects of the treatment according to
the method of claim 1.
29. A method as claimed in claim 28, wherein the sample or another
sample of the subject has been subjected to said analysing before
commencement of the treatment or during diagnosis of the
disease.
30. A method as claimed in claim 28, wherein the treatment
comprises treatment with a drug that treats vessel calcification,
acute myocardial infarction, fracture or stroke.
31. A system, comprising a processor configured for: analysing a
sample of a subject, the sample comprising a material of interest
and at least one adjacent material that is generally adjacent to
the material of interest, the material of interest having a
generally different density from the adjacent material, the
analysing comprising: (i) automatically defining a plurality of
regions of interest within an image of the sample that includes a
junction between the material of interest and the adjacent
material, each of said regions of interest having a width of one or
more pixels or voxels; (ii) determining respective profiles of
density, intensity or attenuation within the regions of interest,
each of the respective profiles crossing the junctions between the
material of interest and the adjacent material; (iii) automatically
identifying the material of interest in the profiles; (iv)
analysing said profiles by automatically defining, for each of said
profiles, a series of sections of said respective profile, each of
said sections comprising at least one S-like portion of the
respective profile, the S-like portion containing the junction
between the material of interest and the adjacent material; and (v)
locating the junction in each of the S-like portions, including:
(a) identifying the greatest difference in values of the respective
profile between an adjacent peak and trough in a segment of the
respective profile that includes at least some of the material of
interest and at least some of the adjacent material; and (b)
locating a point of inflexion in said segment.
32. A system as claimed in claim 31, configured to quantify the
effects of a treatment of a subject, to assess a structure of bone
during a treatment in a subject, to diagnose or monitor a disease
in a subject, or to diagnose or monitor effects of a disease in a
subject.
33. A non-transient computer readable medium containing program
instructions for causing a computer to perform the method of claim
1.
34. A computing device provided with program instructions that,
when executed by the computing device or by a processor of the
computing device, cause the computing device or processor of the
computing device to perform the method of claim 1.
Description
RELATED APPLICATION
[0001] This application is a Continuation of U.S. patent
application Ser. No. 15/191,079, filed 23 Jun. 2016, which is a
Continuation of U.S. patent application Ser. No. 14/720,004, filed
22 May 2015, which is a Continuation of U.S. patent application
Ser. No. 13/395,379, filed 24 May 2012, which is a National Stage
Application of PCT/AU2010/001181, filed 10 Sep. 2010, which claims
benefit of Serial No. 2009904407, filed 11 Sep. 2009 in Australia
and which applications are incorporated herein by reference. To the
extent appropriate, a claim of priority is made to each of the
above disclosed applications.
FIELD OF THE INVENTION
[0002] The present invention relates to an image analysis method
and system for identifying (including automatically if desired) an
object or structure from its surrounding environment, separating or
distinguishing an object or structure from its surrounding
environment, and analysing the object or structure, using--for
example--an image of the object or structure, and is of particular
but by no means exclusive application in the analysis of biological
tissue including the analysis of bone (in vivo or otherwise) and
diagnosis of disease such as bone disease, and the analysis or
characterization of materials and geological samples. In a
particular embodiment, the invention relates to the identification
and distinguishing of bone from surrounding tissue (such as muscle
and fat), the analysis of its structure and the detection of
decayed or abnormal bone (and hence in some case
fracture-vulnerable bone). In another particular embodiment, the
invention relates to the location or identification of a foreign
body embedded within biological tissue (such as muscle).
BACKGROUND OF THE INVENTION
[0003] Currently, bone density is the commonly used tool to predict
fracture risk and evaluate the effects of treatment on bone but it
lacks sensitivity and specificity. More than half of all hip
fractures occur in women determined to be without osteoporosis with
dual energy X-ray absorptiometry (DXA) with a T-score.ltoreq.-2.5.
Furthermore, most women with osteoporosis as defined by DEXA BMD do
not sustain fractures (Delmas P. D., Seeman E. Changes in bone
mineral density explain little of the reduction in vertebral or
nonvertebral fracture risk with anti-resorptive therapy, Bone 34(4)
(2004) p. 599-604). During drug therapy, only a small proportion of
the fracture risk reduction is explained by the increase in bone
density (Ott S. M., When bone mass fails to predict bone failure,
Calcif. Tissue Int. 53(Suppl1) pp. S7-S13). This lack of
sensitivity and specificity partly is partly due to the fact that
fracture risk is not only determined by bone density but also bone
structure (Zebaze R. M. and Seeman S., Measuring Femoral Neck
Strength: Lost in Translation?, IBMS BoneKEy 5 (2008) pp.
336-339).
[0004] Consequently, techniques such as quantitative computed
tomography (CT) and magnetic resonance imaging (MRI) are being used
to produce images that are then analysed to predict fracture risk,
quantify the effects of drugs and determine the effectiveness of
various strategies (such as exercise or diet) for the prevention of
bone loss. However, parameters derived from these imaging
modalities have not proved to be significantly better than bone
density measurement, and cannot presently be used as a substitute
for bone density measurement.
[0005] The problem with these approaches lies not only in image
acquisition but also, and more significantly, in image processing.
One of the most important problems with existing techniques for
assessing bone structure from such images is their reliance on
fixed arbitrary thresholds to identify bone within a region of
interest (ROI) and to segment (i.e. separate) bone into its various
compartments (usually cortical and trabecular bone). Bone
architecture differs from person to person and, for this reason,
the use of fix thresholds has different and unpredictable
consequences on bones from different individuals. Hence, the use of
a fixed threshold can create apparent differences where none exist,
or obscure differences when they exist.
[0006] Furthermore, current methods employ a simplistic model of
bone structure. In particular, these methods treat bone as a
structure made of two distinct compartments (cortical and
trabecular bone). In fact, the so-called `compact` bone and
trabecular (spongy') bone are extremes of a continuum of variation
in porosity from the periosteum to the inner marrow cavity.
[0007] In addition, current methods employ so-called `phantom
calibration` to determine the density associated with each pixel
within an image. This density is then used to identify bone and
assess its structure (hence the term quantitative computed
tomography (QCT)). Existing calibration procedures are specific to
each scanner manufacturer, so bone structure analysis--performed by
software embedded in the scanner--is manufacturer specific.
SUMMARY OF THE INVENTION
[0008] According to a first broad aspect of the invention, there is
provided a method for quantifying the effects of a treatment of a
subject, comprising: [0009] analysing a sample of the subject, the
sample comprising a material of interest and at least one adjacent
material that is generally adjacent to the material of interest,
the material of interest having a generally different density from
the adjacent material, the analysing comprising: [0010] (i)
automatically defining a plurality of regions of interest within an
image of the sample that includes a junction between the material
of interest and the adjacent material, each of said regions of
interest having a width of one or more pixels or voxels; [0011]
(ii) determining respective profiles of density, intensity or
attenuation within the regions of interest, each of the respective
profiles crossing the junctions between the material of interest
and the adjacent material; [0012] (iii) automatically identifying
the material of interest in the profiles; [0013] (iv) analysing
said profiles by automatically defining, for each of said profiles,
a series of sections of said respective profile, each of said
sections comprising at least one S-like portion of the respective
profile, the S-like portion containing the junction between the
material of interest and the adjacent material; and [0014] (v)
locating the junction in each of the S-like portions, including:
(a) identifying the greatest difference in values of the respective
profile between an adjacent peak and trough in a segment of the
respective profile that includes at least some of the material of
interest and at least some of the adjacent material; and (b)
locating a point of inflexion in said segment.
[0015] In an embodiment, the method further comprises applying the
treatment to the subject. This may include performing the analysing
of the sample before, during and/or after the treatment.
[0016] In another embodiment, the treatment comprises treatment
with a drug that treats vessel calcification, acute myocardial
infarction, fracture or stroke.
[0017] In one embodiment, the treatment comprises treatment with a
drug that treats bone and/or modifies bone architecture. The
treatment may comprise an anabolic therapy. The method may comprise
determining any differences in effects of different treatments on
bone architecture, and/or assessing an effect of the treatment on
bone using at least one result of the analysing of the sample of
bone.
[0018] In an embodiment, the method further comprises diagnosing a
disease including the analysing of the sample.
[0019] According to a second broad aspect of the invention, there
is provided a method of assessing a structure of bone during a
treatment in a subject, the method comprising: [0020] treating the
bone; and [0021] analysing a sample of bone of the subject, the
sample comprising a material of interest that includes the
structure and at least one adjacent material that is generally
adjacent to the material of interest, the material of interest
having a generally different density from the adjacent material,
the analysing comprising: [0022] (i) automatically defining a
plurality of regions of interest within an image of the sample that
includes a junction between the material of interest and the
adjacent material, each of said regions of interest having a width
of one or more pixels or voxels; [0023] (ii) determining respective
profiles of density, intensity or attenuation within the regions of
interest, each of the respective profiles crossing the junctions
between the material of interest and the adjacent material; [0024]
(iii) automatically identifying the material of interest in the
profiles; [0025] (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and [0026] (v) locating the junction in each of
the S-like portions, including: (a) identifying the greatest
difference in values of the respective profile between an adjacent
peak and trough in a segment of the respective profile that
includes at least some of the material of interest and at least
some of the adjacent material; and (b) locating a point of
inflexion in said segment.
[0027] In an embodiment, the method comprises assessing an effect
of the treatment using at least one result of the analysing of the
sample of bone. In an embodiment, the method comprises assessing
bone fracture risk of the subject.
[0028] In another embodiment, the method comprises one or more of:
[0029] i) detecting normal bone by analysing the structure of the
material of interest; [0030] ii) detecting diseased bone by
analysing the structure of the material of interest; and [0031]
iii) detecting decayed or abnormal bone by analysing the structure
of the material of interest.
[0032] In a certain embodiment, the method comprises detecting
decayed or abnormal bone by analysing the structure of the material
of interest, including quantifying cortical architecture, wherein a
cortical compartment comprises compact-cortex, outer transitional
zone, inner transitional zone or any combination thereof. In one
example, the method includes: i) quantifying cortical architecture
by determining porosity, cortical thickness, cortical area,
cortical density; or ii) quantifying cortical architecture by
determining porosity, cortical thickness, cortical area, cortical
density and distinguishing cortical porosity attributable to
age-related bone decay from normal porosity.
[0033] In another embodiment, the method comprises detecting
decayed or abnormal bone by analysing the structure of the material
of interest, including quantifying trabecular architecture. In an
example, the method includes quantifying trabecular architecture by
determining trabecular density, trabecular number, trabecular
thickness, trabecular number, trabecular bone volume fraction.
[0034] In an embodiment, the method comprises detecting decayed or
abnormal bone: [0035] i) by analysing the structure of the material
of interest, including combining cortical and trabecular
architecture; [0036] ii) by analysing the structure of the material
of interest, wherein the decayed or abnormal bone comprises both an
abnormal cortex and an abnormal trabecular compartment; or [0037]
iii) by analysing the structure of the material of interest,
including determining any differences in effects of different
treatments on bone architecture.
[0038] According to a third broad aspect of the invention, there is
provided a method of diagnosing or monitoring a disease or the
effects of a disease in a subject, the method comprising: [0039]
analysing a sample of the subject, the sample comprising a material
of interest and at least one adjacent material that is generally
adjacent to the material of interest, the material of interest
having a generally different density from the adjacent material,
the analysing comprising: [0040] (i) automatically defining a
plurality of regions of interest within an image of the sample that
includes a junction between the material of interest and the
adjacent material, each of said regions of interest having a width
of one or more pixels or voxels; [0041] (ii) determining respective
profiles of density, intensity or attenuation within the regions of
interest, each of the respective profiles crossing the junctions
between the material of interest and the adjacent material; [0042]
(iii) automatically identifying the material of interest in the
profiles; [0043] (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and [0044] (v) locating the junction in each of
the S-like portions, including: (a) identifying the greatest
difference in values of the respective profile between an adjacent
peak and trough in a segment of the respective profile that
includes at least some of the material of interest and at least
some of the adjacent material; and (b) locating a point of
inflexion in said segment.
[0045] In an embodiment, the method comprises diagnosing vessel
calcification, acute myocardial infarction, fracture or stroke. The
method may comprise monitoring a disease or the effects of a
disease, wherein the subject is being treated or has been treated
with a drug that treats vessel calcification, acute myocardial
infarction, fracture or stroke. The method may comprise diagnosing
fracture-vulnerable bone.
[0046] In another embodiment, the method comprises one or more of:
[0047] i) detecting normal bone by analysing the structure of the
material of interest; [0048] ii) detecting diseased bone by
analysing the structure of the material of interest; and [0049]
iii) detecting decayed or abnormal bone by analysing the structure
of the material of interest.
[0050] The method may comprise detecting decayed or abnormal bone
by analysing the structure of the material of interest, including
quantifying cortical architecture, wherein a cortical compartment
comprises compact-cortex, outer transitional zone, inner
transitional zone or any combination thereof. In an example, the
method includes: i) quantifying cortical architecture by
determining porosity, cortical thickness, cortical area and
cortical density; or ii) quantifying cortical architecture by
determining porosity, cortical thickness, cortical area and
cortical density and distinguishing cortical porosity attributable
to age-related bone decay from normal porosity.
[0051] In another embodiment, the method comprises detecting
decayed or abnormal bone: [0052] i) by analysing the structure of
the material of interest, including quantifying trabecular
architecture; [0053] ii) by analysing the structure of the material
of interest, including quantifying trabecular architecture by
determining trabecular density, trabecular number, trabecular
thickness, trabecular number and trabecular bone volume fraction;
[0054] iii) by analysing the structure of the material of interest,
including combining cortical and trabecular architecture; or [0055]
iv) by analysing the structure of the material of interest, wherein
the abnormal bone comprises an abnormal cortex and an abnormal
trabecular compartment.
[0056] In a further embodiment, the method comprises assessing an
effect of a treatment of the subject using at least one result of
the analysing of the sample.
[0057] According to a fourth broad aspect of the invention, there
is provided a method of treatment of a disease in a subject, the
method comprising: [0058] quantifying the effects of the treatment
according to the method of the first aspect.
[0059] In an embodiment, the sample or another sample of the
subject has been subjected to said analysing before commencement of
the treatment or during diagnosis of the disease.
[0060] In one embodiment, the treatment comprises treatment with a
drug that treats vessel calcification, acute myocardial infarction,
fracture or stroke.
[0061] According to a fourth broad aspect of the invention, there
is provided a system, comprising a processor configured for: [0062]
analysing a sample of a subject, the sample comprising a material
of interest and at least one adjacent material that is generally
adjacent to the material of interest, the material of interest
having a generally different density from the adjacent material,
the analysing comprising: [0063] (i) automatically defining a
plurality of regions of interest within an image of the sample that
includes a junction between the material of interest and the
adjacent material, each of said regions of interest having a width
of one or more pixels or voxels; [0064] (ii) determining respective
profiles of density, intensity or attenuation within the regions of
interest, each of the respective profiles crossing the junctions
between the material of interest and the adjacent material; [0065]
(iii) automatically identifying the material of interest in the
profiles; [0066] (iv) analysing said profiles by automatically
defining, for each of said profiles, a series of sections of said
respective profile, each of said sections comprising at least one
S-like portion of the respective profile, the S-like portion
containing the junction between the material of interest and the
adjacent material; and [0067] (v) locating the junction in each of
the S-like portions, including: (a) identifying the greatest
difference in values of the respective profile between an adjacent
peak and trough in a segment of the respective profile that
includes at least some of the material of interest and at least
some of the adjacent material; and (b) locating a point of
inflexion in said segment.
[0068] The system may be configured to quantify the effects of a
treatment of a subject, to assess a structure of bone during a
treatment in a subject, to diagnose or monitor a disease in a
subject, or to diagnose or monitor effects of a disease in a
subject.
[0069] The present invention also provides a non-transient computer
readable medium containing program instructions for causing a
computer to perform the method of the first, third or fourth
aspects. The present invention also provides a computing device
provided with program instructions that, when executed by the
computing device or by a processor of the computing device, cause
the computing device or processor of the computing device to
perform the method of the first, third or fourth aspects.
[0070] It should be noted that any of the various features of each
of the above aspects of the invention can be combined as suitable
and desired. It should also be noted that embodiments that refer to
the analysis of bone may also be applied, suitably modified where
necessary, to other samples.
BRIEF DESCRIPTION OF THE DRAWING
[0071] In order that the invention may be more clearly ascertained,
embodiments will now be described, by way of example, with
reference to the accompanying drawing, in which:
[0072] FIG. 1 is a schematic view of a system for detecting
fracture-vulnerable bone according to an embodiment of the present
invention;
[0073] FIG. 2 is a schematic view of the image processor of the
system of FIG. 1;
[0074] FIG. 3 is a schematic view of the memory of the processing
controller of the image processor of the system of FIG. 1;
[0075] FIG. 4 is another, more detailed schematic view of the image
processor of the system of FIG. 1;
[0076] FIG. 5 is a system flow diagram of the system of FIG. 1
(including interaction with a user and customer);
[0077] FIG. 6 is an exemplary image of bone as output by the image
processor of the system of FIG. 1;
[0078] FIG. 7A is a schematic depiction of the determination of the
centre of mass of a bone sample by the image processor of the
system of FIG. 1;
[0079] FIG. 7B is a schematic depiction of the rotation of an image
of a bone sample by the image processor of the system of FIG. 1,
both clockwise and anticlockwise;
[0080] FIG. 8A is an exemplary image derived by the image processor
of the system of FIG. 1 from a DICOM file of a region of interest
(ROI);
[0081] FIG. 8B is the image of FIG. 8A shown in negative for
additional clarity;
[0082] FIG. 9 is a plot of bone strip density (in mgHA/cc) against
line number for the image of FIG. 8A as determined by the image
processor of the system of FIG. 1;
[0083] FIG. 10 is a plot of the first derivative of the density
profile curve of FIG. 9, as determined by the image processor of
the system of FIG. 1;
[0084] FIG. 11 is a plot comparable to that of FIG. 9, illustrating
the determination of the periosteal boundary by the image processor
of the system of FIG. 1;
[0085] FIG. 12 is a plot comparable to that of FIG. 9, illustrating
the determination of the beginning of the compact cortex by the
image processor of the system of FIG. 1;
[0086] FIG. 13 is a plot comparable to that of FIG. 9, illustrating
the determination of the end of the cortical mass by the image
processor of the system of FIG. 1;
[0087] FIG. 14 is a plot of the first order derivative data of FIG.
10, binarized by the image processor of the system of FIG. 1, for
determining the beginning and the end of the cortex;
[0088] FIG. 15 is a plot comparable to that of FIG. 9, illustrating
the determination of the end of the compact cortex by the image
processor of the system of FIG. 1;
[0089] FIG. 16 is a plot comparable to FIG. 9, summarizing the
results of the identification steps performed by the image
processor of the system of FIG. 1;
[0090] FIG. 17 is another exemplary image derived by the image
processor of the system of FIG. 1 from a DICOM file of a region of
interest (upper register) and a copy of this image shown in
negative for additional clarity (lower register), with compartment
boundaries identified by the image processor of the system of FIG.
1 labelled;
[0091] FIG. 18 is a plot of bone strip density (in mgHA/cc) against
line number for the image of FIG. 17 as determined by the image
processor of the system of FIG. 1;
[0092] FIG. 19 is an exemplary image of a bone sample with
`floating cortical` bone identified with the system of FIG. 1;
[0093] FIGS. 20A and 20B illustrate the preparatory steps performed
by the image processor of the system of FIG. 1 for rotating a
region of interest within an image;
[0094] FIG. 20C is a flow diagram of the image recognition and
image analysis rotations performed by the image processor of the
system of FIG. 1;
[0095] FIG. 21 is a BVE scale used by image processor of the system
of FIG. 1 to analyse and characterize the bone during the analysis
rotation;
[0096] FIGS. 22A, 22B and 22C depict the process of image
binarization performed by the image processor of the system of FIG.
1, where FIG. 22A is an exemplary image of a bone sample, FIG. 22B
illustrates the identification of pores within the cortical mass in
one column, and FIG. 22C is a binarized plot of bone versus
non-bone pixels along the line analysed;
[0097] FIG. 23 is a flow diagram illustrating the use of the system
of FIG. 1 to identify an object of known dimensions within an
image, isolate and determine its dimensions;
[0098] FIG. 24 is a schematic depiction of the determination of the
geometric centre of a bone sample by the image processor of the
system of FIG. 1 for the case where the pixel of representative
density is not generally central to the sample;
[0099] FIG. 25 is a flow diagram illustrating the use of the system
of FIG. 1 to identify an object of known dimensions within an
image, isolate the object and determine its dimensions;
[0100] FIG. 26 is a view of a control screen of the image processor
of the system of FIG. 1;
[0101] FIGS. 27A, 27B, 27C, 27D, 27E and 27F are plots, resulting
from the exemplary use of FIG. 25, of thickness, radius, external
perimeter, cross section, internal perimeter and internal surface
area as a function of the ROIs respectively;
[0102] FIG. 28 is a cross-sectional grayscale image of an exemplary
phantom sample for use with the system of FIG. 1, comprising a
plastic tube with drilled holes;
[0103] FIG. 29 is a plot of porosity of several plastic tubes
comparable to that shown in FIG. 28 obtained using the system of
FIG. 1 against known porosity;
[0104] FIG. 30 is a plot of porosity of four phantom samples of
different but known concentrations of hydroxapatite obtained using
the system of FIG. 1 against known porosity;
[0105] FIGS. 31A and 31B are micrographs of bone samples analyzed
with the system of FIG. 1;
[0106] FIGS. 32A and 32B are micrographs of further bone samples
analyzed with the system of FIG. 1;
[0107] FIG. 33 is a micrograph of another bone sample analyzed with
the system of FIG. 1;
[0108] FIG. 34A is a plot of BMD T-scores against porosity for 24
specimens, illustrating the poor correlation between the levels of
porosity measured directly from SEM images and BMD diagnosis
thresholds in the 24 specimens analyzed;
[0109] FIG. 34B is a plot of percentage porosity determined by the
system of FIG. 1 against percentage porosity measured directly from
SEM images;
[0110] FIG. 35 illustrates the positions of points m, n and o
within a density profile curve determined by the system of FIG. 1
when compared with direct observation of the corresponding SEM
image, for a bone sample;
[0111] FIG. 36 illustrates the positions of points m, n and o
within a density profile curve determined by the system of FIG. 1
when compared with direct observation of the corresponding SEM
image, for another bone sample;
[0112] FIG. 37A is a plot of age-related increase in activation
frequency (Ac.F) per year measured directly by histomorphometry as
a function of age;
[0113] FIG. 37B is a plot of age-related increase in porosity
measured using the system of FIG. 1, for comparison with the data
of FIG. 37A;
[0114] FIG. 38A is a greyscale image of a distal radius
reconstructed from high-resolution peripheral quantitative computed
tomography (HRpQCT) data;
[0115] FIG. 38B is an image of the sample depicted in FIG. 38A,
reconstructed by the system of FIG. 1 after extraction of the
periosteal surface;
[0116] FIG. 38C is another view of the sample of FIG. 38A, showing
the segmentation of the same image by the system of FIG. 1 into
four bone compartments;
[0117] FIG. 38D is a plot extracted from FIG. 38C of the outer
transition zone only;
[0118] FIG. 38E is a plot extracted from FIG. 38C of the inner
transition zone only;
[0119] FIGS. 38F and 38G are greyscale versions of, respectively,
FIGS. 38B and 38C;
[0120] FIG. 39A is an image generated by the system of FIG. 1 of a
sample comprising the head of a safety pin embedded in muscle;
[0121] FIG. 39B is an enlargement of the initial region of interest
of the image of FIG. 39A;
[0122] FIG. 40 includes a plot of the density profile curve (upper
register) associated with ROI.sub.1 of the image of FIG. 39A and a
plot of the function .lamda..sub.1 (lower register) associated with
the density profile curve, both generated by the system of FIG.
1;
[0123] FIG. 41A is a cross sectional view of a distal radius
generated by a system of the background art;
[0124] FIG. 41B is a cross sectional view of the same distal radius
generated by the system of FIG. 1;
[0125] FIG. 42A is a cross sectional view of a tibia, with bone and
background delineated by a system of the background art; and
[0126] FIG. 42B is a cross sectional view of the same tibia, with
bone and background delineated by the system of FIG. 1.
DETAILED DESCRIPTION OF THE EMBODIMENTS
[0127] According to an embodiment of the present invention, there
is provided a system for detecting fracture-vulnerable bone, shown
schematically at 10 in FIG. 1.
[0128] System 10 comprises a CT scanner 12, CT control system 14
and an image analyser in the form of image processor 16. Image
processor 16 includes a user interface 18 comprising an input
device and an output device. The input device is in the form of a
keyboard 20 for controlling image processor 16, and the output
device is in the form of a display 22 for displaying images from
the CT scanner 12 before and after processing by image processor
16. CT scanner 12 is configured to perform CT scans of a sample
located within central scanning volume 24 of CT scanner 12 and to
output digitized scan data to CT control system 14; CT control
system 14 is configured to generate image data from the data
received from CT scanner 12, and to output these image data to
image processor 16.
[0129] In this embodiment, the image data comprises a set of
individual image slices or strips through the sample but, in other
embodiments, other sections may be used (such as coronal,
transverse or oblique sections).
[0130] It will be appreciated that system 10 may operate in online
and off-line modes. In the online mode, image processor 16 receives
data directly from CT control system 14 (during or soon after
scanning of the sample). Such an arrangement may be used in a
clinical setting, especially if the detection of any
fracture-vulnerable bone is urgent. In this online mode, the data
is transmitted to image processor 16 via a data link (connected to
respective USBs of the CT control system 14 and image processor 16)
or a communications network (to which CT control system 14 and
image processor 16 are both connected, such as in the form of the
internet); this link or network is shown schematically in FIG. 1 at
26.
[0131] In the off-line mode, image processor 16 receives data
collected earlier by CT scanner 12 and CT control system 14; the
data may be transmitted to image processor 16 via communications
link or network 26, or by any other suitable means (including on
portable computer readable medium, such as CD-ROM, flash card or
the like).
[0132] Image processor 16 includes a processing controller 28 that
is in data communication with input 20 and output 22 and configured
to process image processing instructions in accordance with a
processing procedure (discussed below) and to output processing
outcomes (which may comprise images and/or detection results) to
display 22.
[0133] Referring to FIG. 2, processing controller 28 comprises a
digital processor 30 that processes the processing instructions in
accordance with the processing procedure and--and described
above--outputs processing outcomes to display 22. Keyboard 20 and
display 22 together constitute a user interface 32. Typically, the
processing instructions are stored as program code in a memory 34
of processing controller 28 but can also be hardwired. Herein the
term "processor" is used to refer generically to any device that
can process processing instructions in accordance with the
processing procedure and may comprise: a microprocessor,
microcontroller, programmable logic device or other computational
device, a general purpose computer (e.g. a PC) or a server.
[0134] FIG. 3 shows a block diagram of the main components of
memory 34. Memory 34 includes RAM 36, EPROM 38 and a mass storage
device 40. RAM 36 typically temporarily holds program files for
execution by the processor 30 and related data. EPROM 38 may be a
boot ROM device and/or may contain some system or processing
related code. Mass storage device 40 is typically used for
processing programs.
[0135] FIG. 4 is another schematic view of user interface 32 and
processing controller 28 of system 10 of FIG. 1, with more detail
shown in processing controller 28. Specifically, processor 30 of
processing controller 28 includes a display controller 42 for
controlling display 22, a DICOM file poller 44 for polling DICOM
(`Digital Imaging and Communications in Medicine`) files from a
staging area (typically CT control system 14), a DICOM file
converter 46 for converting DICOM files into images for display, a
file selector 48 controllable to select the files to be processed
by image processor 16, a file digitizer 50 for digitizing the
selected files, a ROI selector 52 for selecting regions of
interest, a bone identifier 54 for identifying bone within a region
of interest, a specimen position adjuster 56, a density profile
analyzer 58 for analyzing the density profile of a sample of bone,
a ROI rejecter 60, an initial result integrator 62, a final result
integrator 64, an index determiner 66, and a result output 68.
Index determiner 66 is adapted to determine various indices
characterizing the bone sample, such that image processor 16 can
transforms image data indicative of the appearance of a bone sample
into indices indicative of various physical properties of the
sample.
[0136] Memory 34 of processing controller 28 includes a ROI result
storage 70 and a file result storage 72.
[0137] The processing instructions that embody the processing
procedure are in the form, as described above, of program code
stored in memory 34, and are adapted to control processing
controller 26 to perform the processing procedure, as shown
schematically in FIG. 5. Referring to FIG. 5, image processor 16
receives (in either online or off-line mode) a DICOM file
(collected--possibly in vivo--of a sample of bone) by polling CT
control system 14, and stores the file in memory 34. The received
file, in this embodiment, is a DICOM file but in other embodiments
may be in a different image format, such as JPEG or TIFF. Image
processor 16 converts the file into an image and outputs the image
to display 22 so that the user can view the image to be processed.
FIG. 6 is an example of such an image 80 of a sample of bone 82
surrounded by soft tissue 84. This allows the user to assess the
adequacy of the sample by verifying that the appropriate file (e.g.
correct bone site) has been uploaded and that the image is of
sufficient quality for the processing to proceed. If these
conditions are satisfied, the user controls image processor 16 to
select the individual files from within the DICOM file (that is,
corresponding to specific image slices through a region of interest
(ROI) that he or she wishes to analyze) for processing. If the user
does not specify the individual files to be processed, image
processor 16 continues with a default selection of the number of
slices. The user then controls image processor 16 to commence the
processing procedure.
[0138] The user may also control system 10 to skip the step of
viewing the image if sufficiently confident that the image is of
high enough quality such that no additional verification is needed.
In these circumstances, system 10 automatically uploads and
analyses the sample without further user intervention, and image
processor 16 initiates analysis without instructions from the
user.
[0139] Initially, image processor 16 digitizes the files using
conventional techniques. Any other desired, conventional
pre-processing may also be performed at this point.
[0140] The first step performed by image processor 16 is to
identify the material of interest (in this example, bone). This
initial identification of the material of interest is done by
localizing a point with a density representative of the density of
the material. In the case of a CT scan image sample containing bone
tissue, a point (B) with the representative density is the pixel
with the highest attenuation (or density). (In other embodiments,
the point with the density representative of density of the
material to be analysed may the point with the lowest attenuation
value, which may be--for example--for the initial identification of
an abscess, as abscesses have low attenuation values.)
[0141] From point B (cf. FIG. 7A), image processor 16 uses radii
r.sub.1, r.sub.2, r.sub.3 and r.sub.4 of the bone at 0.degree., 90,
180.degree. and 270.degree. respectively to determine a geometric
centre (C) of coordinates (x,y), in this embodiment according
to:
( x , y ) = ( ( r 1 + r 3 ) 2 , ( r 2 + r 4 ) 2 ) . ( 1 )
##EQU00001##
[0142] Image processor 16 also determines the centroid (i.e. centre
of mass) G (x.sub.G, y.sub.G) from the initial centre (I) as:
{ x G = i = 1 n A i * x i A i * x i Y G = i = 1 n A i * Y i A i * Y
i , ##EQU00002##
[0143] A.sub.i is the of attenuation (or density) of pixel i and
x.sub.i and y.sub.i its coordinates in the referential of centre
(C). This process is illustrated in FIG. 7A for the example of FIG.
6.
[0144] Depending on the indices to be computed, image processor 16
uses either C or G.
[0145] Image processor 16 next determines the appropriate arm width
(AW) for the analysis. Image processor 16 does not analyse the
entire cross section (or file content) at the same time. Each file
is subdivided into regions of interest (ROI) of specific length
(L.sub.ROI) and width. A ROI may be referred to as an `arm`, so the
width of the ROI is called `arm width` or AW. The AW used by image
processor 16 depends on the indices to be determined and, as is
described below, the size of the material to be analysed. Image
processor 16 employs an appropriate AW in the form of:
AW.ltoreq.r.sub.1*tan(arccos(.DELTA.)), (2)
[0146] where .DELTA. is the estimated curvature of the bone from
the chosen centre (C or G) and r.sub.i is the smallest of the four
radii (viz. r.sub.1, r.sub.2, r.sub.3 and r.sub.4) calculated
above.
[0147] For most indices, image processor 16 commences from centre
C, and .DELTA. is set to a value of 0.98; this AW is designated as
AW.sub.1. For a .DELTA. of 0.98, there is a negligible effect
(0.02) from the curvature of the bone on the analysis. That is, the
bone appears--in effect--essentially locally `flat`, which is why
the arccosine value of 0.98 is chosen. For example, the smallest
radius at the distal radius is .about.8 mm giving image processor
16 an AW of .about.1.62 mm at the distal radius. The smallest
radius at the tibia is 12 mm, for which image processor 16 employs
an AW of .about.2.42 mm. The AW is chosen to be wide enough to
allow pores and resorption cavities to be detected in the density
profile curve (as seen above). Pores and resorption cavities can be
up to 600 microns in diameter in humans, so the AW used by image
processor 16 should be greater than or equal to 1 mm when the
porosity of human bone is being assessed. Notably, as no two bones
have exactly the same dimensions, images of no two bones are
processed in the same fashion by image processor 16, and the
processing of bone structure according to this embodiment is
specific to each bone.
[0148] For other indices, such as strength parameters or external
dimensions, image processor 16 commences from G and uses another
arm width (represented by AW.sub.2). AW.sub.2 is determined by
image processor 16 with a much smaller curvature A (0.999)
according to:
AW.sub.2=r.sub.i*tan(arccos(0.999)). (3)
[0149] AW.sub.2 is smaller than AW.sub.1 because the desire here is
to have the periosteal and enocortical surfaces small enough to be
representative of small arcs, as the immediate aim is to perform an
accurate calculation of the internal and external geometry of the
cross section (rather than detect porosity). At the distal radius
and tibia respectively, AW.sub.2 is between .about.360 microns and
590 microns.
[0150] Image processor 16 can employ a greater default curvature,
but a value of .DELTA. of 0.999 is employed in this embodiment. A
.DELTA. value of 0.999 is used in this example to determine the
precision and the accuracy of the analysis provided by image
processor 16. It should also be noted that .DELTA. is a user
controllable setting in image processor 16. This and other such
settings allow the user to minimize the effect on the analysis of
irregularities of the periosteal rim (i.e. deviation from a portion
of circle). As such it is controllable by the user to minimize any
potential error produced by periosteal rims of ROIs deviating from
portions of circles. The default value of .DELTA.=0.999 is used by
image processor 16 because, with this value, and in view of the
ranges of human bone diameters in vivo (typically less than 5 cm),
AW.sub.2 is 1 mm or less. That is, the periosteal boundary of a
bone is regular at a local level (within, say, a length of 1
mm).
[0151] As image processor 16 employs a calibration, image processor
16 can readily identify muscle tissue and uses muscle density as a
referent. The ratio of muscle to non-muscle tissue is then used as
a form of calibration. This is why, for many calculations,
identifying the muscle and retrieving its density is employed in
the processing performed by image processor 16 (for which reason
image processor 16 employs a different arm width, AW.sub.3, as
described below). Two of the advantages resulting from this are (i)
cumbersome and resource consuming daily calibration for the purpose
measuring a referent density is not needed when scans are to be
analysed using image processor 16, and (ii) image processor 16 can
operate and analyse image files from essentially any CT, not only
quantitative computed tomography. This may obviate the need of
having a QCT. Furthermore, with this approach, image processor 16
can provide existing CT (or MRI) with a bone structure analysis
function. (Existing CT machines, though widely available, lack a
bone analysis function so cannot be used for the diagnosis of
osteoporosis).
[0152] Many indices such as trabecular parameters or indices of
strength or porosity require an a priori detection of the muscle
density, so image processor 16--in these instances--employs a
larger arm width, AW.sub.3:
AW.sub.3=r.sub.x*tan(arccos(0.965)), (4)
[0153] where r.sub.x is the radius selected as the specific angles
outside the bone, with these angles preset based on local anatomy.
Soft tissue outside the bone at these specific angles is
predominantly muscle. Image processor 16 employs, for example,
90.degree. for the radius and .about.65.degree. for the tibia.
AW.sub.3 is chosen to be wide enough to include sufficient muscle
tissue so that its density can be adequately determined, including
for use as a calibration referent.
[0154] After determining the centres and the appropriate arm width,
image processor 16 determines the appropriate angle of rotation
.theta.. Image processor 16 has a default angle of rotation
.theta..sub.1 of 5.degree., which is appropriate for most bones
owing to the expected range of radii. However, this angle may not
be appropriate in some circumstances, so the appropriateness of the
rotation angle is tested before proceeding to the analysis of other
ROIs (as described below). A .theta..sub.1 of 5.degree. gives 72
ROIs per slice analysed.
[0155] For indices using AW.sub.2 as a rotatory arm width, image
processor 16 has a preset angle of rotation .theta..sub.2 of
1.5.degree.. This is small enough for the bone edges (viz.
periosteal surface) to appear locally linear (or at least able to
be approximated by an arc of a circle). This gives a total of 240
ROIs. Image processor 16 is operable to employ a preset angle of
rotation .theta. of 1.degree. (or less) if desired by the user.
[0156] In order to test the appropriateness of the angle of
rotation .theta., image processor 16 verifies whether the following
condition is satisfied:
sin .theta.*r.sub.i<AW, (5)
[0157] where r.sub.i is the smallest of the four radii defined
above.
[0158] If this condition is not satisfied, for analysis using
AW.sub.1, the angle of rotation .theta..sub.1 is deemed unsuitable
and image processor 16 selects the next smaller angle below
5.degree. (such as 4.degree., but in any event suitable for use as
a rotation angle) according to equation 5. This condition ensures
that no fragment of the cross section is left unanalysed. Image
processor 16 informs the user if an angle less than 5.degree. has
been used.
[0159] The same procedure is performed if AW.sub.2 and
.theta..sub.2)(1.5.degree. are used and angles smaller than
1.5.degree. are selected above by image processor 16 as above using
equation 5.
[0160] It should be noted that the angle of rotation (.theta.) and
the AW (as determined by .DELTA.) are parameters set or selected by
the user to control the analysis by image processor 16. Although
they have preset values, they are not fix. The preset values are
chosen for an optimal analysis of bone structure in humans. The
user can thus adjust .theta. and .DELTA. by controlling the
settings controls of image processor 16 for an optimal analysis of
bones, in particular non-humans bones.
[0161] For example, bone diameter can reach up to 10 cm in some
animals, such as horses. To ensure that the length of periosteal
rim remains within a millimetre in length (to minimize the effect
of any potentially curvature irregularity), .DELTA. is set to
0.9999 and .theta..sub.2 readjusted to 0.5.degree.; the larger
cross section of such a bone is subdivided in 480 ROIs of
.about.700 micron length each, allowing an accurate measurement of
the perimeter and surface.
[0162] Image processor 16 then selects an initial region of
interest (ROI) which, in this embodiment, is vertical. This initial
ROI is a transection of the width of image processor 16 AW.sub.1
(or AW.sub.2) starting at C (or G) and of the length going from the
centre and length L, being half of the maximal diagonal of the
entire selection from the centre (C or G), that is, L= {square root
over (a.sup.2+b.sup.2)} (where a and b are defined below). It
should be noted that, in this embodiment, the material (or,
strictly speaking, the image of the material) may be regarded as
within a parallelogram with sides of length a and b.
[0163] When the ROI has been selected, image processor 16 properly
positions the bone within the ROI before proceeding. The process of
positioning involves rotating the bone sample anticlockwise or
clockwise from respectively the first or the last column so that
the beginning of the cortex starts at the same line in the first
and last columns. This ensures that the bone sample is horizontally
positioned. (No rotation is performed by image processor 16 if the
sample is already horizontal.) The suitable positioning of the bone
sample is performed to avoid an angular position of the sample that
would affects the density profile.
[0164] Before actually performing the rotation, image processor 16
identifies the apparent cortical mass of the bone sample. This is
the cortical mass uncorrected for the angular (i.e. non-horizontal)
position of the sample. Within that apparent cortical mass, image
processor 16 determines the beginning of the bone at the first and
at the last columns. This allows the algorithm to recognize if the
sample is perfectly horizontal. As mentioned above, if the sample
is perfectly horizontal no rotation is performed; otherwise, the
appropriate rotation is performed.
[0165] Referring to FIG. 7B (top register), if an anticlockwise
rotation is to be performed, the beginning of the bone at the first
column (i.e. x=0) starts after the beginning of the bone at the
last column (X.sub.a). Image processor 16 detects the beginning of
the cortex at the first column (0, Y.sub.a) and the beginning of
the cortex at the last column (X.sub.a, Y.sub.b) and determines
that Y.sub.a>Y.sub.b; image processor 16 then determines the
angle .delta., between the sample and the horizontal, as:
.delta.=a tan((Y.sub.a-Y.sub.b)/X.sub.a).
[0166] It should be noted that X.sub.a=AW.
[0167] Next, image processor 16 uses .delta. to calculate the
number of pixels (n) by which each column x.sub.i between 0 and
X.sub.a needs to be rotated so that the specimen is horizontal
as:
n=x.sub.i cos(.delta.).
[0168] All columns are then moved anticlockwise by n pixels.
[0169] It should be noted that, because image processor 16 deals
with small values of AW (i.e. relative to the radius of the bone),
viz. less than 3 mm, the angle .delta. is so small that a
tan(.delta.).about..delta. and the number (n) of pixels by which
each column x needs to be rotated clockwise and anticlockwise is
approximated as:
n=x.sub.i*(Y.sub.a-Y.sub.b)/AW
where x.sub.i is the ith column.
[0170] Referring to FIG. 7B (bottom register), if a clockwise
rotation is to be performed, the beginning of the cortex at the
first column (i.e. x=0) starts after the beginning of the cortex at
the last column (X.sub.a). Image processor 16 detects the beginning
of the cortex at the first (0, Y.sub.a) and the beginning of the
cortex at the last column (X.sub.a, Y.sub.b) and determines that
Y.sub.a<Y.sub.b. Next, image processor 16 determines angle
.delta. between the sample and the horizontal as:
.delta.=a tan((Y.sub.a-Y.sub.b)/X.sub.a).
[0171] Using .delta., image processor 16 then calculates the number
of pixels (n) by which each column x.sub.i between 0 and X.sub.a
needs to be rotated so that the specimen is horizontal as:
n=x.sub.i cos(.delta.).
[0172] All columns are moved anticlockwise by n pixels
[0173] Image processor 16 then analyses the density profile strip
by strip. A strip is a rectangular section of the length of the ROI
(i.e. AW, such as AW.sub.1 or AW.sub.2) and of the width of one
voxel (or pixel). (Image processor 16 is also controllable to
modify the characteristics of the strip). FIG. 8A is an exemplary
image derived from a DICOM file from CT scanner 12 of a ROI, with
the individual lines constituting the image labelled from 1 (in
fact outside the bone sample) to 70 (within the bone). FIG. 8B is
the same imagine in negative, for clarity.
[0174] After performing the adequate rotation to adjust the
position of the ROI, image processor 16 then computes a density
profile within the ROI. The type of density profile computed
depends on the purpose of the analysis. Depending on the analysis,
a profile of maximal, minimal or median values per strip may be
computed. However, in most cases, image processor 16 computes a
profile of mean values per strip. (For example, the profile of
maximal values may be computed by image processor 16 when it is
desired to separate two or more bones within an image. Furthermore,
in alternative embodiments it may be desired to control image
processor 16 to calculate, for example, the average square of all
values in defining a profile.) An example of density profile
determination and analyses is described below for the case of mean
values. (The same analysis would be performed if density profile of
maximal values, for example, was computed).
[0175] Image processor 16 determines the mean density DS.sub.i of
each strip i of bone as follows:
D S i = i = 1 m A i / m ( 6 ) ##EQU00003##
[0176] where A.sub.i is the attenuation value of each voxel of the
m voxels within strip i. A curve of the density profile is created
using the mean density profile for each strip. Image processor 16
also identifies the strip with the maximum density (DS.sub.max)
within the entire curve. (At the voxel or pixel level, the words
density and attenuation may be used interchangeably.) FIG. 9 is a
plot of bone strip density (in milligrams of hydroxyapatite per
cubic centimetre or mgHA/cc) against strip or line number, for the
image of FIG. 8A, with the location of maximum strip density
DS.sub.max indicated.
[0177] Then, a first order (or quasi-) derivative of the density
profile curve is calculated by computing the differences in
densities between consecutives strips as follows:
.DELTA.DS.sub.i,i+1=(DS.sub.i-DS.sub.i+1)/DS.sub.i. (7)
[0178] (This definition of derivative, it will be noted, is the
negative of the conventional definition.) From this first order
derivative, image processor 16 analyses the bone sample by
characterizing the sample as comprising a series on `hills` and
`valleys`. Hills correspond to the decreasing parts of the density
profile curve (i.e. where the first order derivative, as defined in
equation 7, is positive) and valleys correspond to the increasing
parts of the density profile curve (where this first derivative is
negative). The derivative has the value 0 at the junctions of the
hills and valleys.
[0179] The derivative of the density profile curve of FIG. 9, as
determined by image processor 16, is plotted in FIG. 10. Image
processor 16 determines the `heights` of the hills as the
difference between the density of the top of the hill and the
density at the bottom of the hill. Image processor 16 determines
the `heights` of the valleys as the difference between the density
at the end of valley and the density at the beginning of the
valley, the opposite of the calculation for hills. That is, hills
are calculated as 1-0 (the density corresponding to 1 before the
series of 0 values minus the density corresponding to the last 0 in
the series), whereas the valleys are 0-1 (the density corresponding
to 0 before the series of 1 values minus the density corresponding
to the last 1 in the series). Using these hills and valleys, image
processor 16 identifies key points in the attenuation profile using
the following method, such as the highest peaks and troughs in the
density profile curve within the cortical mass.
[0180] Image processor 16 then identifies the periosteal boundary j
(i.e. where the bone begins). The density profile curve is
retrieved from the matrix and provided with a new referential of
centre P.sub.1 (cf. FIG. 11). P.sub.1 has a horizontal coordinate
equal to that of DS.sub.max (determined from the derivative of the
density shown in FIG. 10) and a vertical coordinate equal to
minimum density external to DS.sub.max. In this reference frame, a
function .lamda..sub.1(N.sub.i)--being the distance function from
P.sub.1 to the density profile curve--is determined as follows:
.lamda..sub.1(N.sub.i(x.sub.i,y.sub.i))= {square root over
(x.sub.i.sup.2+y.sub.i.sup.2)}.
where N.sub.i refers to strips between the first strip in the
density profile and DS.sub.max (i.e. the strip with the maximum
attenuation). The point j corresponds to the minimum value of
.lamda..sub.1.
[0181] FIG. 11 is a plot comparable to that of FIG. 9, but
illustrating this step.
[0182] Image processor 16 then identifies the beginning of the
cortex by examining the valleys in the first derivative of density
(see FIG. 10) between j and DS.sub.max; the beginning of the cortex
k is the minimal value of the first order derivative in the highest
(or deepest) valley between J and DS.sub.max. This point is
indicated at point k FIG. 10.
[0183] Next, image processor 16 identifies the beginning of the
compact cortex (l), determined as the minimum of the function
.lamda..sub.2:
.lamda..sub.2(N.sub.i(X.sub.i,y.sub.i))= {square root over
(x.sub.i.sup.2+y.sub.i.sup.2)}.
[0184] This step is illustrated in FIG. 12. Identifying l this way
allows the elimination of artefacts such as an eventual local
non-linearity of the periosteum or partial volume effects. In the
identification of l, the density profile curve is reanalysed from a
new referential of centre P.sub.2, with N.sub.i limited to the
strips between k (the beginning of the cortex) and DS.sub.max.
[0185] It will be noted that .lamda..sub.2 is identical in form
with .lamda..sub.1. They differ (as do .lamda..sub.3 and
.lamda..sub.4 discussed below) by the referential point from which
they operate and the portion of the density profile curve to which
they are applied.
[0186] As illustrated in FIG. 13, image processor 16 identifies the
end of the cortical mass (o) (and hence the beginning of trabecular
bone), defined as the minimal value of the function
.lamda..sub.3:
.lamda..sub.3(N.sub.i(X.sub.i,y.sub.i))= {square root over
(x.sub.i.sup.2+y.sub.i.sup.2)}.
[0187] To identify o, the density profile curve is reanalysed from
a new referential of centre P.sub.3, with N.sub.i limited to the
strips between DS.sub.max and z (the end of the curve,
corresponding to the end of the ROI). It will be appreciated by the
skilled person that the order in which these points are determined
is not according to their order in the density profile curve but
rather according to the order in which they are needed for
subsequent steps. This is why o, for example, is determined before
points m and n (see below).
[0188] Next, image processor 16 identifies the end of the
trabecularised cortex n, by identifying the maximal value of the
first derivative in the highest hill between DS.sub.max and o. This
point is indicated at point n in FIG. 10.
[0189] To facilitate the identification of k and n, image processor
16 may `binarizes` the first order derivative of the density (as
shown, for example, in FIG. 10); the result of this binarization is
plotted in FIG. 14. This facilitates the determination of the
greatest of the hills and valleys within the cortical mass (i.e.
between points j and o), which are indicated at 150 and 152
respectively. (The value of the density of each of the valleys and
hills is shown in units of mgHA/cc.)
[0190] Image processor 16 identifies the end of the compact cortex
(m), determined as the minimum of the function .lamda..sub.4 (as
illustrated in FIG. 15):
.lamda..sub.4(N.sub.i(x.sub.i,y.sub.i))= {square root over
(x.sub.i.sup.2+y.sub.i.sup.2)}.
[0191] In identifying m, the density profile curve is reanalysed
from a new referential of centre P.sub.4, with N.sub.i limited to
the strips between DS.sub.max and n (the end the trabecularized
cortex).
[0192] FIG. 16 summarizes the results of these identification steps
performed by image processor 16. Referring to FIG. 16, image
processor 16 has by this point identified in the density profile
(reading from left to right in the figure):
[0193] j: beginning of the bone (i.e. the periosteal boundary);
[0194] k: the beginning of the cortex;
[0195] l: the beginning of the compact cortex;
[0196] m: the end of the compact cortex;
[0197] n: the end of the trabecularised cortex; and
[0198] o: the end of the cortical mass (and hence the beginning of
trabecular bone.
[0199] Once these points are defined, image processor 16 resolves
the bone into four compartments:
[0200] (i) The compact (or `hard`) cortex: between l and m;
[0201] (ii) The trabecularized cortex: between m and n;
[0202] (iii) The cortico-trabecular junction: between n and o;
and
[0203] (iv) Trabecular bone between o and the end of the sample
(z).
[0204] FIGS. 17 and 18 present another example of an image and
density curve determined by image processor 16. FIG. 17 includes
the image (upper register) and--for clarity--negative of the image
(lower register). FIG. 18 is a plot of the determined density
curve, with the Current threshold of separation between cortical
and trabecular bone indicated at 180. Compartment boundaries, as
determined by image processor 16, are indicated for both the image
and the density curve in these figures.
[0205] The analytic segmentation (i.e. separation) of bone has the
advantage of not relying on the use of thresholds. For example, in
existing systems a threshold of 550 to 560 mgHA/cc is used to
separate cortical from trabecular mass (Sharmilla). The approximate
position of this threshold is indicated in FIG. 16 at 160; hence,
after the initial rise in density, existing techniques will deem
all bone internal to the point where the density subsequently drops
below this threshold to constitute trabecular bone. In this
example, the bone mass contained in box 162, though identified as
cortical mass by image processor 16, would thus be identified
incorrectly as trabecular mass by existing methods (even though in
this example essentially no trabeculae are seen, as is apparent
from the absence of data points in the trabecular compartment with
density greater than that of muscle density). The portion of cortex
(indicated at 164) with density above 550 mgHA would be seen by
existing approaches as floating bone in the marrow space (see
`floating cortical` bone 192 in bone image 190 in FIG. 19) after a
threshold has been used to segment the bone into cortical and
trabecular bone; this is an artefact of existing approaches created
by threshold based segmentation of bone into cortical and
trabecular compartments. Furthermore, existing approaches would
also fail to identify the large pore visible in the centre of the
cortex of this example, hence producing an incorrectly low value
for the sample's cortical porosity. Furthermore, existing
approaches would underestimate the decay that has occurred in the
trabecular compartment as they will incorrectly identify cortical
bone as trabecular bone, when trabecular bone has indeed been lost,
leading to an incorrect (and in fact an excessively low) estimate
of trabecular bone loss. Image processor 16 avoids these
pitfalls.
[0206] Box 182 in FIG. 18 illustrates another case where most of
the cortex and its porosity is missed and erroneously labelled as
trabecular bone.
[0207] Image processor 16 then rotates the image anti-clockwise and
performs the same analysis in all the ROIs within the cross section
(or slice). The number n of ROIs is a function of the rotation
angle (.theta. in degrees) as:
n=360/.theta.. (8)
[0208] The method for rotating from ROI.sub.1 to ROI.sub.n is
performed by image processor 16 as follows (and as shown
schematically in FIGS. 20A and 20B), for each slice
[0209] At this point, image processor 16 creates a virtual mirror
matrix (see FIG. 20B, lower register), which has the (x,y)
coordinates of each pixel relative to the centre (G or C) but no
attenuation value. This virtual matrix is the rotation matrix. An
initial region of interest (ROI.sub.1) is then selected from this
rotation matrix. It should be noted that ROI.sub.1 is selected so
that its length is half that of the longest diagonal, so that no
point with an attenuation value in the initial matrix is lost
during the rotation process. Image processor 16 then retrieves the
attenuation value of each point (or pixel) in the virtual matrix
corresponding to ROI.sub.1 from the initial matrix.
[0210] Image processor 16 then rotates each point A(x,y) in
ROI.sub.1 of this virtual matrix to A'(x',y') according to:
x ' = x cos .theta. + y cos ( .theta. + .pi. 2 ) and ##EQU00004## y
' = x sin .theta. + y sin ( .theta. + .pi. 2 ) . ##EQU00004.2##
[0211] Once rotated, the attenuation value of pixel A' is retrieved
from the initial matrix. Image processor 16 rotates the whole
initial ROI (ROI.sub.1) at the same time. Image processor 16
rotates every points in the initial ROI.sub.1 and repopulates a
matrix similar in width and length to ROI.sub.1 until a full
360.degree. has been performed. For example, for an angle .theta.
of 5.degree., each cross section is analysed as 72 matrices or 72
ROIs. For an angle .theta. of 1.5.degree., the entire cross section
is analysed as 240 matrices (or ROIs); for an angle .theta. of
0.1.degree., the entire cross section is analysed as 3600 matrices
(or ROIs).
[0212] Image processor 16 then ignores any strips with points with
no corresponding attenuation values in the initial matrix owing to
a particularly long ROI.sub.i.
[0213] After selecting the appropriate AW and angle (.theta.),
image processor 16 analyses the image by performing a variable
number of rotations around the entire cross section. For most
indices (such as porosities, mineralization and trabecular
parameters), a minimum of two rotations is employed in the analysis
of the image.
[0214] FIG. 20C is a schematic flow diagram 200 of the operation of
image recognition and image analysis stations of image processor
16. At step 202, the image is inputted by image processor 16 into
the image recognition station.
[0215] At step 204, it is determined whether a calibrant is
required. For the measurement of parameters (e.g. cortical
thickness, external dimensions and internal dimensions) that do not
require a calibrant such as the dimensions of the material,
recognition rotations are not performed. Thus, if a calibrant is
not required, processing by image processor 16 passes to step 208
(see below); otherwise, processing continues at step 206.
[0216] At step 206, one or more `recognition` rotations--as many as
required--are performed, to identify characteristics or features
within the image that will be used in subsequent rotations for a
detailed quantitative analysis. These features include one or more
representative densities of surrounding soft tissue and one or more
representative densities of the material of interest (in this
example, bone) in the sample.
[0217] Thus, in the case of bone within an image, the objective of
the recognition rotation(s) at step 206 is to (i) identify the
surrounding background including the muscle tissue and retrieve a
representative density of, say, the muscle tissue that will be used
as calibrant in subsequent analysis rotations, and (ii) identify a
representative density of material to be analysed, which will also
be used subsequently in the quantitative analysis of the sample.
Processing then passes to step 208.
[0218] At step 208, image processor 16 inputs the image into the
image analysis station, to analyse and characterize the material of
interest. At step 210, image processor 16 employs the image
analysis station to perform one or more analysis rotations (as many
as required) to analyze the sample. During the analysis rotation(s)
at step 210, image processor 16 uses the characteristics retrieved
during the image recognition rotations (see step 206), in this
example bone, to quantify the parameters as described above.
[0219] At step 212, image processor 16 exits the image analysis
[0220] station and begins results synthesis.
[0221] The analysis performed by image processor 16 exploits
differences in density from one strip to the next within a
compartment, so any compartment (viz. compact cortex,
trabecularized cortex, cortico-trabecular junction and trabecular
bone) limited to one strip is unsuitable for this analysis as, in
this example, porosity, decay, etc., cannot be determined by image
process 16. In addition, a compartment limited to one strip is
unanalysable because the confounding effect of partial volume
effects (PVE) on this strip cannot be assessed. Hence, image
processor 16 does not analyse compartments limited to one
strip.
[0222] Image processor 16 views a cortical mass of four strips (or
indeed fewer) as unsuitable for analysis: the three intracortical
compartments will have only one strip (which are hence each
unanalysable) and the last or fourth strip may be tainted by
PVE.
[0223] However, regardless of the number of strips within the
cortical mass, image process 16 approximates the thickness of
cortex (CThcm) as follows:
CTHcm = ( n - 1 ) * re + ( DS n * re ) DS max , ##EQU00005##
[0224] where n is the number of strips within the cortical mass, re
is the resolution, DS.sub.n and DS.sub.max the densities of the
last strip n and the density of the strip with the maximal density
in the adjacent suitable ROI.sub.i-1 or ROI.sub.i+1. (If both ROIs
ROI.sub.i-1 and ROI.sub.i+1 are suitable, the DS.sub.max is
retrieved from the ROI with highest DS.sub.max). It is assumed
here, that last strip (n) is tainted by PVE.
[0225] In addition, internal and external radii are perimeters can
be calculated as the beginning of the bone can be identified and
the thickness of the cortex can be approximated.
[0226] Using this approach, image processor 16 selects and analyses
n ROIs of the same width as ROI.sub.1. Once all the ROIs within the
slice have been analysed, image processor 16, then merges the n
ROIs to reconstruct the analysed slice that has the characteristics
required for the recognition or the diagnosis.
[0227] It should be noted that, during the rotation step, ROIs of
interest overlap; that is, a pixel (i) may appear in more than one
ROI. For example, at the distal tibia, for an angle .theta. of
0.1.degree., the cross section is reconstructed with 3600 ROIs and
the same point may appear in as many as 90 ROIs. When this is so,
image processor 16 takes into consideration (such as by tracing
back a point in the ROIs to the same point in the original matrix)
the contribution of overlapping points to each specific parameter,
as it is important that the contribution of each pixel should
counted only once. Image processor 16 counts and stores the number
of times a point appears in a specific zone (viz. compact cortex,
trabecularized cortex or trabecular bone) in different ROIs. During
the merging process, image processor 16 first determines the
frequency distribution of the appearance of the point in different
zones, then attributes the point to the zone in which it appears
most frequently (though other criteria may be used to determine the
correct or most appropriate zone).
[0228] After reconstructing the first slice, image processor 16
advances to the next slice until all slices or all the preselected
slices have been analysed. When all the slices have been analysed
image processor 16 determines various parameters and outputs the
relevant ones for the individual sample.
[0229] Image processor 16 allows a very broad investigation of bone
and its surrounding tissues (typically comprising muscle and fat).
The application of image processor 16 is not limited to indices
relevant to the detection of osteoporosis. Indices produced
relevant to many other metabolic diseases and infiltrative bone
diseases such bone metastasis. Image processor 16 is controllable
by the user to choose the indices the user desires, and to perform
the chosen analysis promptly.
[0230] Cortical areas and thicknesses are determined by image
processor 16 from the geometric centre C, with the arm width
AW.sub.1 and the rotatory angle .theta..sub.1. These indices, when
output, are displayed in units of millimetres or centimetres (for
thicknesses) and mm.sup.2 or cm.sup.2 (for surfaces).
[0231] The average of the cortical thicknesses of the compact
cortex (CTh.sub.cc) is analogous to the mean CTh provided by
existing QCT techniques, but is determined differently. According
to this embodiment, image processor 16 determines a true average
according to:
i = 1 n CTh i / n , ( 9 ) ##EQU00006##
[0232] where n is the number of ROIs analysed for the individual,
that is, the mean of all CThs values (cf. current QCT techniques,
which determines an estimate of the mean computed as Cortical
area/Cross sectional perimeter).
[0233] The value determined by image processor 16 is the average of
the 72 regions (or more) analysed in n slices. Image processor 16
allows the user to output the standard deviation (SD) and standard
error of the mean (SEM) and the range for this value. Hence, the
outputs are the average CTh, SD, SEM.
[0234] The median value of cortical thicknesses of the compact
cortex is not determined in existing approaches, but the present
inventors have ascertained that median CTh provides a better
prediction of bone strength as estimated with finite element
analysis than the average Cth, so image processor 16 determines and
outputs median CTh. Image processor 16 determines the median CTh as
the median of all 72 (or more) ROIs.times.n slices measured values
of CTh. Image processor 16 can also determines the range (viz.
minimal and maximal values) of such median values.
[0235] The minimum value of cortical thicknesses of the compact
cortex is not provided by existing QCT techniques, but the present
inventors have concluded that it is likely to be useful in clinical
and research settings. There is evidence that the minimal value of
CTh within the cross section shows the greatest decrease with age
and therefore is likely to provide a better indication of the
fracture risk that the mean CTh (which is understandable from a
biomechanical perspective, just as a chain breaks at its weakness
link). Image processor 16 can also determine a mean, SD, SEM and
range of values of the minimum value of cortical thicknesses of the
compact cortex, and identify the anatomical location of the compact
cortex (e.g. anterior, posterior, etc) (see the results summarized
below).
[0236] The lower the area of the compact cortex (or `cortical
area`), the greater the risk of fracture. Cortical area (or CoA) is
determined by image processor 16 as the number of pixels within the
compact cortex multiplied by the area of the pixels, hence:
CoA=n(re.sup.2), (10)
[0237] where n is the number of pixels within the compact cortex
and re is the resolution of CT scanner 12 (or other scanner used to
image the bone, such as a MRI scanner). This can also be determined
from:
CoA = i = 1 i = 360 / .theta. ( ( .theta. / 360 ) * ( ( .pi. * r i
2 ) - ( .pi. * ( r - CTh i ) 2 ) ) ) , ( 11 ) ##EQU00007##
[0238] where i is the ROI, n is the number of ROIs examined,
CTh.sub.i is the cortical thickness at ROI.sub.i and rotation angle
.theta. and r.sub.i is the radius of the bone at the ROI.sub.i
(determined as described below). Image processor 16 effectively
divides the cross section into small arcs of radii of curvature
r.sub.i, allowing it to accurately determine areas of irregular
cross section (such as bones). The mean, SD, SEM and range of CoA
can then be determined for all the cross sections examined.
[0239] Image processor 16 is also configured to determine
CTh.sub.tc the average of the cortical thicknesses (CTh) of the
trabecularized cortex. (CTh.sub.tc is not determined in existing
techniques.) CTh.sub.tc is indicative of bone loss and fracture
risk independently of the average thickness of the compact cortex.
A thicker trabecularized cortex indicates that more decay has
occurred and is more fragile (even if the residual compact cortex
apparently remains as thick as that of peers). This value is
determined by image processor 16 using equation 10 or 11, but
applied within the compact cortex.
[0240] Image processor 16 is configured to determine the median
value of the cortical thicknesses (CTh) of the trabecularized
cortex, which provides additional information on the status of the
cortex above and beyond the average value, and to determine the
maximum of the cortical thicknesses (CTh) of the transitional (or
equivalently the trabecularized) cortex, which shows where the
cortex is weakest.
[0241] Image processor 16 is also configured to determine a novel
parameter, the area of the trabecularized cortex, which is the area
occupied by the pixels within the trabecularized cortex. It should
be interpreted differently from the area of the compact cortex (of
which it is, in a sense, the complement): the greater the area of
trabecularized cortex, the more the magnitude of cortical decay
that has occurred and hence the greater the fragility. This value
is determined by image processor 16 again using equation 10 or 11,
but applied within the trabecularized cortex.
[0242] Image processor 16 is also configured to determine the novel
index CTh.sub.cm: the average of the cortical thicknesses (CTh) of
the cortical mass. Existing approaches identify the cortical mass
but only in the compact cortex. The thicker the cortical mass, the
better it is. This is determined using equation 9, but applied to
the compact mass.
[0243] Image processor 16 is configured to determine two other
novel indices: the median value of cortical thicknesses of the
cortical mass, and the minimal value of cortical thicknesses of the
cortical mass and their respective anatomical locations (see the
results summarized below). Existing approaches do not identify the
cortical mass, being limited to the compact cortex. The greater the
median value of the cortical mass, the stronger the bone, and the
smaller the minimal value of the cortical mass, the weaker the
bone.
[0244] Image processor 16 is also configured to determine the area
of the cortical mass, as the lower the cortical area, the greater
the risk of fracture. Cortical area is determined by image
processor 16 according using equation 10 or 11, but applied within
the compact mass.
[0245] Image processor 16 is also configured to determine the
average of the cortical thicknesses (CTh) of the cortico-trabecular
junction (represented as CTh.sub.ctj), the median value of cortical
thicknesses of the cortico-trabecular junction, and the minimum
value of cortical thicknesses of the cortico-trabecular junction.
These indices are of value as, while in healthy bones there is a
clear differentiation between the cortical and the trabecular
compartment so the cortico-trabecular junction is very small, with
aging and decay the cortico-trabecular junction increases so is
blurred.
[0246] Image processor 16 is similarly configured to determine the
area of the cortico-trabecular junction, using equation 10 or 11,
but applied within the cortico-trabecular junction.
[0247] Image processor 16 can determine the respective radii
r.sub.i for each of the 72 or more ROIs defined as the length
between the start of the cortex and (C), as follows:
r.sub.i=(z-1)re, (12)
[0248] where r.sub.i is the radius at ROI.sub.i and re is the
resolution of CT scanner 12.
[0249] Image processor 16 can also determine porosity and other
related indices, such as the percentage compact cortex (PCC), which
is the proportion of the cortex that has a compact or solid
appearance, determined from the ratio of the thicknesses of the
compact cortex to that of the cortical mass calculated as
follows:
PCC=100CTh.sub.cc/CTh.sub.cm. (13)
[0250] A small value of PCC indicates that the cortex has
diminished owing, for example, to ageing. (Certain disease
processes such as hyperparathyroidism may also be indicated).
[0251] Similarly, image processor 16 is configured to determine
percentage trabecularized cortex (PTC), the proportion of the
cortex that has a trabecular-like appearance and is calculated as
follows:
PTC=100CTh.sub.tc/CTh.sub.cm. (14)
[0252] A high PTC means that the cortex has diminished during
ageing or disease, such that a greater proportion of the cortex
resembles trabecular bone.
[0253] Image processor 16 is configured to determine percentage
cortico-trabecular junction (PCTJ), the proportion of the cortex
that is a transitional zone, as follows:
PCTJ=100CTh.sub.ctj/CTh.sub.cm. (15)
[0254] Generally, the higher the porosity, the weaker the bone, so
image processor 16 is configured to determine several porosity
indices.
[0255] For example, image processor 16 is configured to determine
the apparent porosity (which is sometimes referred to simply as the
`porosity`, but which--as the skilled person will
appreciate--should be distinguished from the absolute porosity) of
the compact cortex (aPoCC). If there has been no increase in
porosity, all the strips would have a similar density. Apparent
porosity is therefore calculated as the ratio of (i) the area of a
rectangle defined by DS.sub.max and the thickness of the cortex and
(ii) the area under the density curve within the specific
compartment of the cortex (cf. FIG. 16 for the boundaries l, m, n
and o of the respective compartments). The apparent porosity of the
compact cortex is thus determined according to:
aPoCC ( % ) = 100 ( i = 1 n ( ( DS i + DS i + 1 ) / 2 ) ) / ( n - 1
) DS max , ( 16 ) ##EQU00008##
[0256] where n is the number of strips within the cortex in the
ROI. This is an approximation of the more general formulation:
aPoCC ( % ) = .intg. x = 1 n DS ( x ) . x DS max * ( n - 1 ) ,
##EQU00009##
[0257] which reflects the fact that the analysis is based on the
ratio of two areas and is not necessarily a point by point
analysis. The more general formulation may be more suitable in
other embodiments.
[0258] Image processor 16 is also configured to determine the
porosity (again, actually apparent porosity) of the trabecularized
cortex (aPoTC) (according to equation 16 but limited to the
trabecularized cortex), of the cortical mass (aPoCM) (according to
equation 16 but for the entire cortex), and of the trabecular
compartment (aPoTB) (according to equation 16 but for the
trabecular compartment).
[0259] Porosities (that is, `absolute` rather than apparent
porosities) are also determined by image processor 16. To quantify
porosities, image processor 16 identifies the representative
densities of the surrounding muscle tissue (Max.sub.md) and a
representative density A(x) of the material of the material of
interest (i.e., bone in this instance)
[0260] To identify Max.sub.md, a minimum of rotation is needed.
From the first rotation, using AW.sub.3, image processor 16
identifies the surrounding muscle tissue and retrieves the maximal
density of the muscle tissue (Max.sub.md). To do so image processor
16 determines the frequency distribution of the attenuation values
of all the voxels corresponding to the surrounding muscle tissue.
Max.sub.md correspond to the attenuation corresponding to the 95th
percentile in the frequent distribution curve. The maximal values
may be tainted by noise from, for example, Compton scattering. The
value at the 95th percentile is large enough to be representative
of the Max.sub.md but essentially untainted by noise.
[0261] The same process is used by image processor 16 to identify
(B) the value representative of pure bone. There is experimental
evidence to suggest that the maximal density of each ROI fulfils
these criteria. The inter-haversian distance exceeds 200 .mu.m, so
four haversian canals define a square area exceeding 4,0000
.mu.m.sup.2 where no pores are present. Hence, at a spatial
resolution of 82 .mu.m, up to six pixels that are not tainted by
partial volume effects will be placed between haversian canals
within the cortex. It can therefore be assumed that the pixel with
the highest density within a given ROI in is indeed of these pixels
and has the density of purely mineralised bone. Hence, for every
imaging modality of 200 .mu.m resolution or below, the attenuation
of (B) is considered as made of purely mineralized bone, untainted
by PVE. Image processor 16, then used the attenuation of (B) is
then used in conjunction with Max.sub.md to convert the attenuation
values of all the voxels or pixels with an image into an scale
ranging from -1 to +1. This transformation is a new concept that
allows the attenuation of each voxel to convert into the estimated
amount of bone within each voxel. The transformed unit of each
voxel is referred to as bone volume equivalent (BVE). FIG. 21 is a
BVE scale used by image processor 16 to analyse and characterize
the bone during the analysis rotation. Fat has a BVE of -1, water
(corresponding, in bone, to cells, blood vessels and interstitial
tissue, hence similar to the surrounding muscle tissue or water)
has a BVE of 0, and pure bone as a BVE greater than or equal to
0.9. Any voxel with a BVE lower than 0.9 is partly (if BVE>0 but
<0.9) or fully (if BVE.ltoreq.0) made of void.
[0262] Image processor 16 calculates the BVE of a voxel (X) of
coordinates (i,j) within the ROI as:
BVE ( X ) = ( A ( X ) - A ( Max md ) ) ( A ( B ) - A ( Max md ) ) ,
##EQU00010##
[0263] where A(X), A(B) and A (Max.sub.md) are respectively, the
attenuation values of the voxels X, B and Max.sub.md.
[0264] The concept of BVE allows image processor 16 to provide a
particularly detailed characterization of the bone, including
estimating porosity within the image due to pores of size below the
nominal in plane resolution of the imaging modality.
[0265] When the image is captured by an imaging modality with an
in-plane resolution of 200 .mu.m or below, voxels with a BVE equal
to or greater than 0.9 are voxels with purely mineralized bone (as
differences in density or attenuation due differences in tissue
mineralization do not exceed 10%). Hence, any voxel with a BVE
between 0 and 0.9 contains a proportion of void. If the voxel is
located within the cortical mass, this voxel contains a pore of
size smaller than the resolution of modality (cf. FIG. 21).
[0266] Voxels with a BVE of 0 are voids of size similar to or
greater than the resolution of the imaging modality.
[0267] Image processor 16 calculates the pore volume equivalent
(PoVE) of each voxel from the BVEs as follows:
[0268] if BVE(X).ltoreq.0, then PoVE(X)=1, and
[0269] if BVE(X)>0, then PoVE(X)=1-BVE(X).
[0270] The void content of each voxel within the image is then
estimated from the PoVE.
[0271] (iii) Image processor 16 then calculates the total PoVE
(TPoVE) within the ROI as:
TPoVE = x = 1 n PoVE ( x ) . ##EQU00011##
[0272] (iv) Finally, image processor 16 calculates the porosity
(Po) within the ROI as:
Po = 100 * TPoVE n , ##EQU00012##
[0273] where n is the number of voxels within the ROI.
[0274] By this procedure, the contribution of each voxel to
porosity is accounted for by image processor 16 proportionally to
estimated proportion of purely mineralised bone within the voxel
using the voxel's density.
[0275] It should be noted the range of pores to be assessed can be
selected by the user with user interface 18 of image processor 16.
For example, pores with PoVE=0 (i.e. voxels that appear to be empty
to within the resolution of system 10), or pores with PoVE=0.5
(i.e. pores of half the size of the resolution of system 10).
[0276] The user can choose the quantify the porosity due to pores
of any size within any of the three compartments of bone (i.e.
compact cortex, trabecularized cortex and trabecular bone).
[0277] A method of quantification of pores of the size of the
resolution of machine is as follows. Image processor 16 converts
each pixel (or voxel) in each compartment into a binary value: bone
`1` and non-bone `0`, with bone defined as pixel i, such that:
A.sub.i>Max.sub.md, (18)
[0278] and non-bone defined as pixel j, such that:
A.sub.j.ltoreq.Max.sub.md. (19)
[0279] The result of this process is shown schematically in FIGS.
22A, 22B and 22C. FIG. 22A is an example of a ROI with a pore
visible within the cortex. After identification of compartments,
pores are defined as voxels with density inferior to Max.sub.md. An
example of the identification of pores within the cortical mass in
one column is shown in FIG. 22B. FIG. 22C is a binarized plot of
bone versus non-bone pixels along the line of analysis.
Binarization and identification of pores occurs after the cortical
mass has been identified in the analysis of the entire ROI.
[0280] Existing techniques for determining porosity use a
dichotomous approach to distinguish pore from bone tissue, based on
a fixed arbitrary threshold. Porosity so determined is highly
resolution dependent and influenced by partial volume effects
(PVE), so image processor 16 regards it merely as `apparent
porosity`. Voxels with variable proportion of bone and soft tissues
at the edge of the pore-bone (i.e. tainted by PVE), such as point i
(which is tainted by PVE) shown in FIG. 22B are erroneously
excluded from the calculation of porosity. The same errors applies
when current algorithms computes porosity after using a threshold
to separate bone from soft tissues. This is why image processor 16
advantageously computes porosity as described (cf. equations 16 and
17) where the density of every voxel is accounted for.
[0281] It should be noted that the general principle embodied by
image processor 16 is the determination of porosity by employing
the surrounding environment and a representative density of
material of interest in its pure form as referents. This enables
image processor 16 readily to provide a non-destructive measurement
of porosity of any material, not only bone. The porosity of
materials (for example, rocks, metals and plastics) may be
determined using a vacancy or `emptiness` as a referent. This is
illustrated, by way of example, in FIG. 23. FIG. 23 is a schematic
view of a porosity measurement apparatus 220 for the
non-destructive in vitro measurement of the porosity of materials,
according to an embodiment of the present invention, for use in
conjunction with system 10 of FIG. 1.
[0282] Apparatus 220 includes a sealable tube 222 with an internal
station 224 (for supporting a material or sample to be analysed
226). Tube 222 has a valve 228, and apparatus 220 includes a vacuum
pump 229 attached to valve 228. Vacuum pump 229 is provided for
evacuating tube 222 once material or sample 226 is on station 224
and tube 222 has be sealed.
[0283] Thus, after mounting material or sample 226 on station 224
and sealing tube 222, vacuum pump 229 is switched on and any
floating material around material or sample 226 is evacuated from
tube 222 by vacuum pump 229. Thus, an at least partial vacuum is
created around material or sample 226, which serves as the vacancy
or `emptiness` (of low or zero density) exploited according to the
present embodiment. Tube 222 can then be detached from vacuum pump
229 and located within central scanning volume 24 of CT scanner 12,
scanned, and the resulting data forwarded to image processor 16 for
processing as described above.
[0284] Establishing this `emptiness` (essentially a region of low
density) around material or sample 226 allows image processor 16 to
use the surrounding environment (in this case the so-called
`emptiness`) to determine porosity within material or sample 226. A
material floating or otherwise around material or sample 226 would
otherwise be identified by image processor 16, which would use it
as a referent, leading to inaccurate results. (Emptiness' is used
in quotes to acknowledge that a truly zero pressure cannot be
achieved, and in any event is not necessary; the removal of this
`floating debris` is sufficient.)
[0285] After impurities (e.g. this floating material) have been
removed, the material of interest is imaged (CT) as described
above. Images are automatically retrieved and analysed by image
processor 16 to determine the absolute porosity using the
surrounding environment (in this case `emptiness`) as referent as
described above. The use of a referent (emptiness in the case of
porosity) avoids the reliance on absolute CT-values.
[0286] In one embodiment, the invention provides an integrated
porosity measurement system comprising system 10, sealable tube 222
and vacuum pump 229.
[0287] It should be noted that in the case of porosity
quantification in a living person, pores are non-bone spaces within
the bone envelope, rather this so-called `emptiness`. These
non-bone spaces are occupied by tissues of density equal or
inferior to that of the muscle. This why assessment of porosity in
bones, in living individuals (in in vivo settings), uses muscle as
referent.
[0288] Non-destructive measurement of porosity has many potential
applications including material characterization, detection of
internal flaws and cracks within the material.
[0289] It should be noted that image processor 16 can output this
porosity, such as for comparative purposes. In this case, the pore
area (PoA) in each compartment is:
PoA = i = 0 n PoVE . ( 20 ) ##EQU00013##
[0290] Image processor 16 then determines porosity in each
compartment. For example the porosity of the compact cortex (PoCC)
is:
PoCC = 100 * PoACC AreaCC , ( 21 ) ##EQU00014##
[0291] where PoCC is the porosity of the compact cortex, PoACC is
the pore area of the compact cortex and AreaCC is the area of the
compact cortex.
[0292] Image processor 16 is configured to determine the decay of
the compact cortex (CD1). This index reflects a gradual increase in
porosity from the inner to the outer cortex as observed during
ageing. Pores are uniformly distributed within the cortex of young
healthy individuals, so--for them--that the line of best fit within
the cortex is substantially horizontal. With ageing, porosity
erodes the inner more than the outer aspect of the cortex, so the
line of best fit lines makes an increasingly greater angle with the
horizontal. A normal CD1 is less than 3.5%.
[0293] Image processor 16 determines CD1 by determining the line of
best fit within the compact cortex using conventional least mean
squares formulae, then determining the angle (.omega..sub.1, in
radians) between the line of best fit and the horizontal line
passing through DS.sub.max, and then determining CD1 as:
CD1(%)=100(.omega..sub.1/1.5707). (22)
[0294] Image processor 16 determines the decay of the
trabecularized cortex (CD2), which is analogous to CD1. CD2 is
calculated by determining the line of best fit within the end of
the compact cortex (m) and the beginning of the trabecular bone
using conventional least mean squares formulae, then determining
the angle (.omega..sub.2, in radians) between the line of best fit
and the vertical line passing through m. CD2 is then determined
as:
CD2(%)=100(.omega..sub.2/1.5707).
[0295] Image processor 16 can also determine cortical fragility
(CF), which is indicative of--and hence allows the assessment
of--the contribution of porosity due to high remodelling to
cortical fragility. CF is determined as follows:
CF (%)=(AbsPo.sub.cc*CD.sub.1+AbsPo.sub.tc*CD2)/10000. (23)
[0296] Image processor 16 can also determine relative density, that
is, the attenuation of a respective compartment (compact cortex,
trabecularized cortex, etc) relative to that of the surrounding
muscle. Relative density thus reflects both the mineralization and
the porosity and is expressed as a percentage. Relative densities
are independent of DS.sub.max.
[0297] The average density (or attenuation) DS.sub.m of muscle
strips is:
DS m = ( i = 1 k - 1 DS i / n ) . ( 24 ) ##EQU00015##
[0298] The relative density of each compartment is then determined
in a similar manner. For example, the relative density of the
compact cortex is:
RD.sub.cc=100DS.sub.cc/DS.sub.m, (25)
[0299] where DS.sub.cc is the density of the compact cortex. A
relative density thus indicates on average how dense the bone is
compared to the surrounding tissues. The less mineralised is the
bone, the less porous, and the lower is its relative density.
[0300] Likewise, the relative density of trabecularized cortex
(RD.sub.tc), that is, the ratio of the density of strips of the
trabecularized cortex (DS.sub.tb) to DS.sub.m, is:
RD.sub.tc=100DS.sub.tc/DS.sub.m.
[0301] The relative density of the cortical mass (RD.sub.cm), that
is, the ratio of the density of strips of the cortical mass
(DS.sub.cm) to DS.sub.m, is:
RD.sub.cm=100DS.sub.cm/DS.sub.m.
[0302] The relative density of the trabecular mass (RD.sub.tb),
that is, the ratio of the density of strips of the cortical mass
(DS.sub.tb) to DS.sub.m, is:
RD.sub.tb=100DS.sub.tb/DS.sub.m.
[0303] The cortico-trabecular differentiation index (CDTI) is the
ratio of the area under the density profile curve in the trabecular
bone compartments to the area under the density profile curve in
the cortical mass compartments, as is determined by image processor
16 from:
CTDI = 100 i = o + 1 z ( ( DS i + DS i + 1 ) / 2 ) i = k o ( ( DS i
+ DS i + 1 ) / 2 ) , ( 26 ) ##EQU00016##
[0304] where k, o and z are respectively the strips corresponding
to the beginning of the compact cortex, the end of the cortical
mass and the end of the entire sample. It should be noted that, in
largely cortical sites (such the midshaft of the femur or the
subtrochanteric region) in young adults, the cortex is clearly
distinct from the trabecular bone. In these sites, bone loss
transforms the cortex into a inner trabecular-like structure and an
therefore the differentiation between the cortex and the trabecular
bone diminishes, thus decreases CTDI.
[0305] Image processor 16 is also configured to determine various
fat indices. These indices reflect the extent of fat transformation
of the bone marrow. The more the bone has been replaced by fat, the
more bone has been lost and the greater the fragility of the bone.
The attenuation in a compartment within the bone that is lower than
the surrounding muscle (of representative density Max.sub.md)
reflects a loss of bone within that compartment and its replacement
by fat. Fat voxels therefore constitute strips within the
trabecular bone that have density less than the muscle density
Max.sub.md.
[0306] Image processor 16 is configured to determine fat proportion
(FatP) (the number of fat voxels divided by the total number of
examined voxels) and the fat burden (FatB) (the average density of
fat voxels (FS) compared to the representative density of muscle,
Max.sub.md). The average density (DS.sub.f) of fat strips is:
DS f = i = 1 n ( D i ) / n , ( 27 ) ##EQU00017##
[0307] where D.sub.i is the density of voxel i and
D.sub.i.ltoreq.Max.sub.md. The fat burden (FatB) is therefore:
FB (%)=100(DS.sub.f/Max.sub.md). (28)
[0308] Image processor 16 is also configured to determine various
mineralization indices. These indices reflect the degree of
mineralization and the heterogeneity in the distribution of the
mineralization. It is assumed that the difference in the
attenuation between the bone pixel with the maximal density and the
muscle (Max.sub.md), which is the contrast between these two
tissues, reflects the level of mineralization of the bone. This
index is not calculated using strips but bands. Image processor 16
identifies pixels or voxels that are almost certainly entirely made
of purely mineralised bone and hence, minimally tainted by partial
volume effects and porosity. These are voxels with a BVE equal to
or greater than 0.9, as discussed above. The mineralization (ML)
level is determined after examining all the cross sections within
the image, not at ROI.
[0309] Image processor 16 determines ML as:
ML ( % ) = 100 * 1 n i = 1 n ( ( A ( i ) - ( A ( Max md ) ) A ( B )
, ( 29 ) ##EQU00018##
[0310] where A(i) is the attenuation of the pixel or voxel (i) such
that BVE(i)>0.9.
[0311] When desired, image processor 16 can determine and output
the absolute densities of the compact cortex, trabecularized
cortex, cortico-trabecular junction and trabecular bone in mgHA/cc
(HA=hydroxyapatite) or mg/cc of mineralised bone, as required by
the user.
[0312] Image processor 16 converts the attenuation of each voxel
within a given compartment of the bone to mgHA/cc with the
equation:
D (mgHA/cc)=(695.808*(A.sub.i/Max.sub.md)))-455.27, (30)
where A.sub.i is the attenuation of voxel i. This equation was
established experimentally by scanning phantoms.
[0313] Image processor 16 converts the attenuation A of each voxel
within a given compartment into its equivalent of mineralised bone
in g/cc, firstly by identifying the voxels associated with muscle
and then using the muscle density as a referent to determine the
density of bone. The need for daily calibration is avoided owing to
the use of the subject's own muscle as a referent. The density (D)
of other voxels in g/cc is then determined by image processor 16
from:
D (g/cc)=0.495*(A/Max.sub.md). (31)
[0314] Purely mineralised bone has a density of 2.14 g/cc.
[0315] Image processor 16 outputs the absolute densities (in g/cc
and mgHA/cc) of cortical and trabecular bone. This is for
comparison only, as absolute densities are more commonly
required.
[0316] Image processor 16 determines strength indices by using the
centroid G, the rotatory arm width AW.sub.2 and the rotation angle
.theta..sub.2. Non-bone tissue within the subperiosteal envelope is
removed using Max.sub.md+1 as discussed above. Image processor 16
then determines strength indices as follows.
[0317] Image processor 16 determines buckling ratios (BR) as
follows for each ROI:
BR.sub.i=CThcc.sub.i/r.sub.i, (32)
where CThcc.sub.i and r.sub.i are respectively the thickness of the
compact cortex and the radius at the ROI.sub.i.
[0318] Image processor 16 determines the second moment area (or
moment of inertia I) for all angles .theta. within the slice as
follows:
I x = .intg. i = 1 n y i 2 ( re 2 ) , ( 33 a ) ##EQU00019##
[0319] where I.sub.x is the second moment of area about the axis x,
re is the resolution of the pixel, and y.sub.i is the perpendicular
distance from the axis x to the element pixel.
[0320] Image processor 16 calculates the second moment area in all
directions around the cross section. However, before calculating
the moment of inertia, image processor 16 excludes any non-bone
pixels from the image. Again, no threshold is used to identify
non-bone pixels. Non-bone pixels are pixels of attenuation below
that of the beginning of the bone (i.e. muscle or below). This
exclusion of non-bone pixels (i.e. pores) allows image processor 16
to calculate I as free as possible from the effects of
porosity.
[0321] The mass adjusted second moment area is calculated in a
similar fashion to the moment of inertia, but with the density or
attenuation of each pixel taken into account. The mass adjusted
second moment area is thus determined as follows for each ROI:
I xma = .intg. i = 1 n ( ( A i ) / A max ) y i 2 ( re 2 ) , ( 33 b
) ##EQU00020##
[0322] where I.sub.xma is the mass adjusted second moment area
about the axis x, A.sub.i is the density (or attenuation) of each
bone pixel and A.sub.max is the maximal attenuation of all pixels
within the slice (i.e. the pixel almost certainly made of purely
mineralised bone). Image processor 16 determines the section
modulus similarly to the determination of I in all directions,
again excluding non-bone pixels, as follows:
z x = .intg. i = 1 n y i 2 re 2 ( .intg. i = 1 n y i ) / n . ( 34 )
##EQU00021##
[0323] Image processor 16 determines the mass adjusted section
modulus similarly to the determination of I in all directions,
again excluding non-bone pixels, as follows:
z ma = .intg. i = 1 n ( ( A i ) / A max ) y i 2 re 2 ( .intg. i = 1
n y i ) / n . ##EQU00022##
[0324] Image processor 16 determines the product moment of area,
excluding non-bone tissue and not taking the angle into account
(which is irrelevant in this case), as follows:
Ixy = .intg. i = 1 n x i y i re , ( 35 ) ##EQU00023##
[0325] where x and y are the coordinates in a frame of reference
with centre G as defined by bone pixels. The product moment of area
is a determinant of bending stress in an asymmetric cross section.
This is the case for most bone, which are asymmetric. Unlike the
second moments of area, the product moment may give both negative
and positive values. From this image processor 16 determines the
maximum and minimum mass adjusted second moments area, as well as
their orientation within the cross section.
[0326] Image processor 16 determines the mass adjusted product
moment area after separating bone from non-bone pixels, and without
taking the angle into account as it is again not needed. Image
processor 16 determines the mass adjusted product moment area as
follows:
S i x = .intg. i = 1 n ( A i ) / A max ) x i y i re , ( 36 )
##EQU00024##
[0327] where A.sub.max and A.sub.i are respectively the maximal
attenuation and the attenuation of pixel i. x and y are the
coordinates of each pixel i in a frame of reference with centre G.
Image processor 16 then determines the maximum and minimum mass
adjusted second moment area, as well as their orientation within
the cross section.
[0328] Image processor 16 is configured to determine the polar
moment of inertia J, which predicts the ability of the cross
section to resist torsion, as follows:
J = .intg. i = 1 i = n ( x 2 + y 2 ) re 2 . ( 37 ) ##EQU00025##
[0329] Image processor 16 determines the polar moment of inertia
for every angle, and also determines its maximum and minimum
values.
[0330] Image processor 16 is configured to determine the mass
adjusted polar moment of inertia J.sub.max, as follows:
J max = .intg. i = 1 n ( A i ) / A max ( x 2 + y 2 ) re 2 . ( 38 )
##EQU00026##
[0331] Image processor 16 is also configured to determine the total
cross sectional area (TCSA) and perimeter. For these indices, image
processor 16 uses the centroid G, the rotatory arm width AW.sub.2
and the rotation angle .theta..sub.2. Non-bone tissue within the
subperiosteal envelope is removed using Max.sub.md as referent as
described above.
[0332] For each ROI.sub.1, image processor 16 computes the
L.sub.ROIi (length of ROI.sub.i) according to:
L.sub.ROI.sub.i=.theta..sub.2*(2ri*.pi.)/360. (39)
[0333] After rotation by 360.degree., image processor 16 determines
the total perimeter P of the cross section from:
P = i = 1 360 / .theta. 2 ( .theta. ( 2 r i .pi. ) / 360 ) . ( 40 )
##EQU00027##
[0334] Image processor 16 then determines the total cross section
area (or periosteal area) TCSA as:
T C S A = i = 1 360 / .theta. 2 ( .pi. r i 2 ) .theta. 2 / 360. (
41 ) ##EQU00028##
[0335] Image processor 16 then determines the endocortical
perimeter ECP as:
ECP = i = 1 360 / .theta. 2 ( .theta. 2 ( 2 ( r i - CThcc i . )
.pi. ) / 360 ) , ( 42 ) ##EQU00029##
[0336] where CThcc.sub.i is the thickness of the compact cortex in
ROI.sub.i
[0337] Image processor 16 then determines the endocortical area ECA
as:
ECA = i = 1 360 / .theta. 2 ( .pi. ( ( r i - CThcc i . ) 2 ) )
.theta. 2 / 360. ( 43 ) ##EQU00030##
[0338] Image processor 16 then computes the ratio ECA/TCSA, which
is an indicator of fragility.
[0339] Image processor 16 treats each ROI as an infinitesimal
portion of a circle of radius r.sub.i, allowing the determination
of the areas and perimeter of structure of complex shape such as a
bone cross section. This is why AW.sub.2 is very small.
[0340] All non-bone tissue is removed from the trabecular
compartment using Max.sub.md as referent so that, for any ROI, from
strip o to z, only trabeculae are included. As discussed above,
image processor 16 employs a rotatory arm width AW.sub.1 and a
rotatory angle .theta..sub.1. (In fact the choice of AW here is of
no great importance; any rotatory AW will suffice. As discussed
above, image processor 16 has three AWs; AW.sub.1 and AW.sub.2 are
rotatory in the sense that they can move from one ROI to another
one, while AW.sub.3 is static and only serves to assess muscle
density at specific locations.)
[0341] After removal of non-bone tissue as explained above, all
non-bone tissue pixels are given the value 0 and bone tissue pixels
the value of 1. The width (w.sub.i) for the n non-bone voids for
each column i is determined as:
w.sub.i=n*re, (44)
[0342] where re is by the resolution of the scanner. Trabecular
size (Tr.S.sub.i) for column i is calculated as:
Tr . S i = ( ( z - o ) re ) - w i m i , ( 45 ) ##EQU00031##
[0343] where m.sub.i is the number of trabeculae. A trabecular bone
is regarded as a series of consecutive pixels coded as 1. The
presence of a `0` indicates the end of the trabecular.
[0344] Image processor 16 repeats this process for each line j,
where the trabecular size for each line j is:
Tr . S j = AW - w j m j , ( 46 ) ##EQU00032##
[0345] and where m.sub.j is the number of trabeculae in the j line.
It should be noted that a column is a series of pixels in the
direction of the density profile and a line is perpendicular to the
column.
[0346] Finally, image processor 16 determines the Tr.S for the ROI
from:
Tr . S i = ( i = 1 n Tr . S i + j = 1 n Tr . S j ) / ( m i + m j )
. ( 47 ) ##EQU00033##
[0347] Image processor 16 determines and outputs this trabecular
size, as this reflects both the shortening and the thinning of
trabeculae. Existing approaches employ a definition of trabecular
thickness that does not account for trabecular decay due to
shortening.
[0348] Image processor 16 determines trabecular separation in each
line (i) and each column (j) from:
Tr.Sep=(w.sub.i/n.sub.i)+(w.sub.j/n.sub.3), (48)
[0349] where n.sub.i and n.sub.j are respectively the numbers of
voids in i columns and j lines. A void is a series of uninterrupted
0s. The presence of a 1 signals the end of the void.
[0350] Image processor 16 can output various results graphically to
display 20, including: [0351] The cortical thickness distribution
around the cross section; and [0352] The cortical thickness along
the bone at any given angle (e.g. anteriorly).
Summary of Results
[0353] Image processor 16 is operable in various modes:
[0354] (i) Default Mode: The results are summarized for each slice,
then for all the slices as average, minimal, median, standard
deviation for indices other partial perimeters and partial
surfaces. Partial perimeters and surfaces are summed for to obtain
the perimeter or the surface for each slice and average to obtain
the average perimeter or surface for the sample. In this default
mode, all the parameters are available.
[0355] (ii) Flexible Analysis Mode: This mode of analysis allows
the user to control image processor 16 to determine some parameters
at specific anatomical locations (anterior, antero-lateral, medial,
etc). Image processor 16 provides this flexibility in the analysis
because a given parameter may affect bone strength and fracture
risk differently depending on the anatomical region. For example at
the femoral neck, an increase in porosity on the superior aspect
may be related to fracture more than a similar increase in porosity
on the inferior aspect.
[0356] Image processor 16 numbers or indexes regions of interest
within a slice, so each ROI corresponds to a specific anatomical
location. For example, at the distal radius, ROI1 correspond to the
lateral aspect of the radius as image processor 16 rotates
anticlockwise, for angle of 5 degrees for example, ROI.sub.36
corresponds to the anterior aspect, and ROI.sub.54 to the medial
aspect.
[0357] ROIs at specific anatomical locations can then be combined
to provide specific indices at that location within the slice or
(cross section). The number of ROIs (N) to be combined in each
anatomical region is estimated to be:
N=(360/(m*.theta.)),
[0358] where .theta. is the rotation angle, and m is the number of
anatomical regions the slice is to be divided into. For example,
the user can control image processor 16 to divide the slice into 2
halves (e.g. medial and lateral) by selecting a value for m of 2
starting at ROI.sub.1, or m=2 starting at ROI.sub.i such as
i=360/(4*.theta.). It should be noted that the minimal value for m
is 2 (as a value for m of 1 would mean that the entire cross
section is analysed, which switches image processor 16 to Default
Mode), whereas the maximal value of m is 360/.theta..
[0359] It should that, in Default Mode, only parameters for which
local or regional values are meaningful are output. These
parameters includes thicknesses, porosities, mineralization,
trabecular indices are indices of strength moment of inertias.
Parameters that are meaningful for the entire cross section, such
perimeters and areas, are not output in Flexible Mode.
[0360] In Default Mode, the average value for the parameter at a
specified anatomical location is determined for all the slices
within the file.
Interpretation of Results
[0361] System 10 is able to assess the structure of bone and its
changes during various interventions (such as exercise, nutrition
variation or treatment) on bone. This has many applications. For
example, it is suggested that treatment by drugs such as strontium
ranelate, parathyroid hormone (PTH) treatment or odanacatib (in
clinical trial phase 3) increase the diameter of the bone and,
therefore, its bending strength. This can be examined by accurately
quantifying bone dimensions without using threshold as done using
system 10. A threshold based assessment of bone structure performed
according to existing methods may not be optimal for this purpose
as, for example, identification of bone using thresholds may result
in erroneous removal new and poorly mineralised on the periosteal
surface leading to the erroneous conclusion the bone has not been
formed, hence periosteal apposition has not occurred when indeed it
has occurred. This is likely to occur when the effect of PTH
treatment on bone are being assessed. PTH is anabolic and
stimulates the formation new and poorly mineralised bone. This can
also occur when assessing the effects of exercise on bone as
exercise especially during growth stimulates the formation of new
and poorly mineralised bone.
[0362] Furthermore, strontium ranelate is absorbed by the bone
tissue and has a higher atomic number that calcium. Hence,
strontium treatment increases the density or attenuation of bone
voxels. An increase in the density of voxels adjacent to edge of
the bone (periosteum or marrow space) may affect edge detection,
giving the erroneous impression that this drug has changed the
dimensions of the bone. This occurs when a density based fixed
threshold used (according to existing approaches) to examine bone
structure before and after treatment with this agent gives an
ambiguous result regarding its effects on bone structure. The
conversion of all attenuations into a single BVE scale ranging from
-1 to 1 and independent of the absolute attenuation of values of
the bone, and the fact that thresholds are not used by image
processor 16, minimizes the risk that the different effects of the
drugs on attenuation or density of voxels will be confused with
differential effects on bone structure (i.e. dimensions such as
cortical thickness) and porosity. (With existing processing
methods, conclusions about differing effects of two treatments on
bone structure and porosity may be made merely because one
treatment has changed the material attenuation more than another
(Rizzoli R, Laroche M, Krieg M A, Frieling I, Thomas T, Delmas P,
Felsenberg D., Strontium ranelate and alendronate have differing
effects on distal tibia bone microstructure in women with
osteoporosis, Rheumatol. Int. 30(10) (2010) pp. 1341-8).)
Furthermore, this avoids the changes in the density of attenuation
values of voxels resulting in erroneous assessment of bone
architectural parameters such as porosity. For example, pores are
defined as voxels with attenuation below a given value; a drug like
strontium ranelate significantly increase the attenuation value of
bone tissue, so a change in bone tissue attenuation may be confused
with a decrease in porosity. Such erroneous assessment is minimized
or eliminated in image processing performed by image processor
16.
[0363] System 10 is also able to determine the level of decay of
bone and to detect fracture-vulnerable bone with a fair degree of
confidence, based on the following considerations.
[0364] (i) A normal cortex
[0365] The absolute values of CTh and cortical area within entire
cortical mass in any of its sub-compartments are not sufficient for
determining whether the cortex is normal. That is, values higher
than average (for the relevant subject group) does not mean that
decay has not occurred. Similarly, the absolutely values of
porosity and cortical area are not sufficient to determine whether
a cortex is normal, as a cortex may be small owning to low growth
and not owing to bone loss. Similarly a high-porosity is not
necessarily the result of increased intracortical remodelling; it
may also be growth-related. Hence, absolute values are merely
indicative. This is a fundamental shift from current approaches
which view a smaller cortex or area compared to peers as
abnormal.
[0366] Instead, image processor 16 determines the architecture of
the cortex and uses the result to determine whether the cortex is
normal or not. A normal cortex has: [0367] a low PTC (<20% based
on current evidence), suggesting little or no trabecularization,
and [0368] a high PCC (>60% based on current evidence).
[0369] The normality of the cortex is further supported by one or
more of the indices output by image processor 16 that is: [0370]
Little if any decay (CD1 normal) [0371] Normal RDcc [0372] Normal
RDcm [0373] A normal cortical thickness or cortical area [0374]
PoCC and PoCM are low [0375] A low index of fragility
[0376] (ii) A decayed cortex
[0377] A cortex that is of similar or greater thickness than that
of its peers may still be decayed (notwithstanding the
classification used in existing approaches, which regard cortex of
similar or greater thickness than that of its peers as normal).
Some individuals are born with thicker than average cortices so,
even after loosing bone, they may still have a cortex of normal
thickness. Such individuals are thus misdiagnosed by existing
approaches as having normal bones. However, image processor 16 does
not necessarily characterize a thick cortex as normal. Rather,
image processor 16 characterizes cortex (as normal or abnormal)
according to its architecture. An abnormal cortex has a high PTC.
In addition this cortex can be thin with a low PCC with a high
porosity and low relative density, and hence high values of CD1
and/or CD2 and a high fragility.
[0378] (iii) A normal trabecular bone has a normal relative
density, low porosity, a low fat proportion and fat burden.
[0379] (iv) An abnormal trabecular bone will have a low relative
density, high porosity and high fat proportion and fat burden all
these suggesting that trabecular bone has been lost and replaced by
fat.
[0380] (v) An abnormal bone will have an abnormal cortex, an
abnormal trabecular compartment, or both an abnormal cortex and an
abnormal trabecular compartment.
[0381] Image processor 16 also employs porosity indices of each
sub-compartment to make such determinations.
[0382] (i) Image processor 16 determines porosity from the area
under the density profile curve, permitting such determinations in
vivo and ensuring that the contribution of each pixel in each
sub-compartment is accounted for. This may be compared to existing
approaches in which porosity is determined from the number of pores
within a specific bone area, which is limited by the resolution of
the image (bearing in mind that most pores are not visible in such
images) and by the threshold approach (also bearing in mind that
most pixels classified as bone and rejected in the assessment of
porosity may be tainted by PVE, that is, partially on the bone and
partially on a pore). The difference between the area under the
curve and the area defined by DS.sub.max is strongly correlated
with the porosity measured by direct histomorphometry. It has been
found that the increase in porosity with ageing does not occur
uniformly within a specimen, so--according to this embodiment--the
strip with the lowest porosity is used by image processor 16 as a
marker for use in assessing the increase in porosity.
[0383] (ii) Image processor 16 can determine porosity in the
trabecularized cortex and the cortico-trabecular junction. Existing
approaches measure porosity only in the compact cortex, but cannot
accurately assess porosity in the trabecularized cortex and the
junctional zone.
[0384] Image processor 16 can determine an index of the degree of
cortical porosity attributable to age-related bone (cortical decay)
within the compact cortex and compact+trabecularized cortex. When
the trabecularized cortex is made of cortical ruins, decay is
present. Young healthy bones have cortices of varying degree of
porosity but no decay. Old bones with evidence of bone loss also
have bones with varying degree of porosity but high amount of
decay.
[0385] Cortical fragility is the product of cortical porosity and
cortical decay. It is expected that greater porosity and greater
decay will result in greater fragility.
[0386] The proportion of the cortex that is compact (PCC), that is,
the ratio of the thicknesses of the compact/(cortical mass)
expressed as a percentage.
[0387] Thicknesses and areas of each sub-compartment:
[0388] The cortical-trabecular differentiation index (CTDI) permits
differentiation between the cortex and the trabecularization. With
ageing and bone loss, the de-differentiation of between the
trabecular compartment and the cortical complex (i.e. cortex plus
transition zone) occurs. This is a marker of decay and
fragility.
[0389] It should also be noted that image processor 16 outputs most
indices in the form of percentages, which are readily interpreted
by the user. For example, it should be immediately apparent that a
high percentage porosity or a low percentage density is
deleterious. (Existing methods output absolute numbers, which are
generally interpretable only by expert users; for example,
densities are expressed in mgHA/cc so the user must know what
constitutes a low, acceptable and high range in order to interpret
such results.)
[0390] In the above, it was assumed that point B is generally in
the centre of the material under analysis. If not, image processor
16 can be controlled to identify a more central point F (see FIG.
24) within the material, such that:
{right arrow over (BF)}={right arrow over (BI.sub.1)}+{right arrow
over (BI.sub.2)}+{right arrow over (BI.sub.3)}+{right arrow over
(BI.sub.4)}
where I.sub.1, I.sub.2, I.sub.3 and I.sub.4 are points further away
from B in the top left, top right, lower left and lower right
quadrants with respect to B, at the edges or beginnings of the bone
(determined as described above). Image processor 16 then determines
C, the centre of mass of the sample, from F in the same way that it
determines C from point B described above.
Examples
[0391] FIG. 25 illustrates the use of system 10 to identify an
object of known dimensions within an image, isolate and determine
its dimensions with high precision and accuracy. This experiment
was used to test the accuracy and the precision of image processor
16 in identifying and determining the dimensions of a structure
within an image.
[0392] FIG. 26 is a view of a control screen of image processor 16.
Various settings can be selected, as can the region (e.g. bone) to
be analyzed (see the second row).
[0393] Referring to FIG. 25, a plastic tube 230 of known uniform
dimensions was scanned 232 with scanner 12 of system 10. The DICOM
file was retrieved from scanner 12; the resulting image file 234
(indicative of an image 236 of the tube) was generated and
transmitted to image processor 16 for analysis.
[0394] The image 236 was analyzed by image processor 16 using
.theta..sub.2, generating 240 ROIs, and AW.sub.2. Image processor
16 identified G and generated the thickness, radius, external
perimeter, cross section (i.e. total area), internal perimeter, and
internal surface area as a function of the ROIs (from ROI.sub.1 to
ROI.sub.240), which are presented in FIGS. 27A, 27B, 27C, 27D, 27E
and 27F respectively.
[0395] This test serves several purposes. Firstly, the precision of
the analysis by image processor 16 can be tested. As the object 230
was known to have a uniform structure, any difference in structural
parameters from ROI to ROI is largely attributable to noise in
image processing. It can be seen that the dimensions (thickness,
radius, etc) of the object were constant (hence almost perfectly
flat lines as plotted) from ROI to ROI. This confirms that the
analytical procedures used by image processor 16 to rotate and
retrieve elements of a structure from ROI to ROI are sound, with
little noise. (If elements of the structure differed from ROI to
ROI, this would indicate that the analysis of a structure from ROI
to ROI by image processor 16 were noisy and poorly reproducible).
The precision of the analysis by image processor 16 was then
calculated as expressed as coefficient (CV) of variation, that is,
the mean of the 240 ROIs divided by the standard deviation of the
240 values. The CV for thickness, radius, external perimeter,
internal perimeter, external surface and internal surface
measurements were respectively 3%, 1.7%, 1.7%, 2.1%, 3.4% and 4.2%,
confirming the high reproducibility of the analysis.
[0396] Secondly, as the dimensions the plastic tube were well
known, the accuracy of image processor 16 in determining the
dimensions of a structure could be evaluated. The size of the
plastic tube was measured directly, and the corresponding
parameters were determined with system 10. The results are
tabulated and compared below:
TABLE-US-00001 Determined by Parameter Measured directly system 10
Agreement thickness of plastic 3 cm 2.9 cm 96.7% radius 19.75 cm
19.99 cm 98.7% external perimeter 124.09 cm 125.6 cm 98.8% internal
perimeter 105.2 cm 106.8 cm 98.5% external surface 881.4 cm.sup.2
908.6 cm.sup.2 97% area internal surface area 1225.4 cm.sup.2
1256.1 cm.sup.2 97.5%
[0397] Altogether, this shows that image processor 16 can
automatically detect an object within an image and determine its
dimensions with good precision and accuracy.
[0398] In addition, image processor 16 was able to automatically
identify the plastic tube and robustly determine its dimensions
without using thresholds (as employed in existing techniques).
Unlikely image processor 16, existing techniques require an
external input to perform such analysis; that is, an operator must
reset the density threshold for detection so that the object can be
identified within the image. This requires knowing a priori the
density of the object and/or manually adjusting the image until the
image of the object appears separated from the surrounding
environment. Image processor 16 is adapted to identify a structure
within an image and separate it essentially regardless of its
density.
[0399] In another example, a variety of numbers of holes were
drilled in a series of plastic tubes (such as that shown in cross
section in grayscale image 280 of FIG. 28) to produce different
levels of porosity. These tubes were then imaged using system 10,
images were collected and porosity was quantified using image
processor 16. FIG. 29 is a plot of the porosities determined by
system 10 against the true porosities (i.e. those determined
experimentally from the known characteristics and numbers of holes
in each tube). A high correlation was found, with an r.sup.2 value
of 0.98.
[0400] In a further example, four cylindrical phantoms of different
but known concentrations of hydroxapatite were imaged using CT
scanner. The images were collected with system 10 and analyzed
using image processor 16 to determine the level of mineralization.
FIG. 30 is a plot of the mineralization level (%) from image
processor 16 against the known hydroxyapatite concentration
(mgcm.sup.-1): high correlation between the two is again
evident.
[0401] Twenty-four bone specimens from cadavers were studied using
scanning electron microscopy (SEM) and QCT, and their bone density
was measured in in vivo like conditions (after submersion in saline
and an adequate internal rotation of 15.degree.). Porosity and
other indices were determined directly from micrographs and
estimates after DICOM file from CT scanner 12 using image processor
16. Porosity profile curves (analogous of density profile) were
generated after porosity was measured directly from SEM images.
Points identified by image processor 16 from the density profile
curves were compared with the image scans and SEM micrographs.
[0402] Examples of these results are shown in FIGS. 31A, 31B, 32A,
32B, 33, 34A, 34B, 35 and 36.
[0403] Micrographs from the subtrochanteric region, bone mineral
density (BMD) in these specimens were measured at the trochanteric
region and converted into T-scores using the Geelong reference
range which is the reference recommended by the Australian and New
Zealand Bone and Mineral Society.
[0404] FIG. 31A is a micrograph from image processor 16 of bone
from a 72 year old woman with a modestly elevated intracortical
porosity as can be seeing on the micrograph. The small increase in
porosity is recognized by image processor 16, which quantifies the
porosity as 4.4%. FIG. 31B is a micrograph from image processor 16
of bone from a 90 year old woman with a markedly elevated
intracortical porosity. This markedly increased in intracortical
porosity is also captured by image processor 16 which quantifies
the porosity as 35.1%. The difference in decay between the two
samples is clearly visible, FIG. 31A showing more decay that FIG.
31B. This was confirmed by porosity results determined with image
processor 16 (whereas according to a conventional bone density
tests, these two bones have similar decay, as bone density tests do
not discern the difference in fragility between these two
bones).
[0405] FIG. 32A is a micrograph from image processor 16 of bone
from a 29 years old woman with minimal porosity; image processor 16
quantified the porosity of this sample as 0.9%. FIG. 32B is a
micrograph from image processor 16 of bone from a 72 year old woman
with normal cortex and low porosity; this low porosity is well
captured by image processor 16, which quantified the porosity of
the sample of FIG. 32B as 0.1%. However, according to the BMD test
(the currently used method of diagnosis), the 29 year old woman is
normal (with a T-score of 1.92) whereas the 72 year old woman is
osteopenic (with a T-score of -1.14), suggesting bone fragility in
the older woman when there is none. Image processor 16 thus has the
potential to avoid such cases of misdiagnosis by adequately
detecting decay when it is present.
[0406] FIG. 33 is a micrograph from image processor 16 of bone from
a 67 year old women with a significantly decayed cortex (as can be
seen in the micrograph). Image processor 16 captured this decay and
quantified porosity as 14.7%. However, the woman was classified
according to the BMD test as normal (with a T-score of -0.2). Many
people with a normal bone density test sustain a fracture at some
point; it is likely that this woman would be one of them. It is
expected that image processor 16 analysis would identify these
people who fracture despite having a normal BMD test.
[0407] FIG. 34A is a plot of BMD T-scores against porosity
(expressed as percentage of bone area), and illustrates the poor
correlation between the levels of porosity measured directly from
SEM images and BMD diagnosis thresholds in the 24 specimens
analyzed. This shows that BMD test poorly captured cortical decay.
About 50% ( 7/13) of the specimens with high porosity (>20% of
bone area) consistent with having lost bone, had normal BMD
(T-score>-1). Only 50% ( 2/4) of specimens with low aBMD
(T-score<-2.5 SD) have high porosity (>20%); the remaining
had low aBMD despite having normal porosity. Of the 10 specimens
with T-scores less than -1 (i.e. osteopenic or osteoporotic), only
5/10 (50%) had low porosity and the remaining 50% had high
porosity.
[0408] FIG. 34B is a plot of percentage porosity determined by
system 10 against percentage porosity measured directly from SEM
images. As is apparent from this plot, porosity as quantified by
image processor 16 correlated well with directly measured
porosity.
[0409] FIGS. 35 and 36 illustrate the positions of points m, n and
o within a density profile curve when compared with direct
observation of SEM images for two respective samples. The ability
of image processor 16 to identified points within the density
profile curve as compared to images generated in vivo from DICOM
files and analysed by image processor 16 was presented above
(detailed description of the invention).
[0410] In these two in vitro samples, it can be seen that the
points and hence the dimensions (in particular thicknesses) within
the samples are clearly identified. The decay of the compact cortex
(CD1) and trabecularised cortex (CD2) are also indicated. Referring
to FIG. 35, sample 1 is from a 29 year old woman with a
histologically normal cortex. This was confirmed by image processor
16, which found a CD1 of 1.4% and a CD2 of 4.2%; this corresponds
to a homogeneous cortex followed by a clear transition from the
cortex to marrow space. The cortex is almost entirely compact
(PPCC=87.5%), the PTC and PCTJ are limited to one strip each, and
therefore inexistent as discussed in the section on unsuitability
of the analysis. The normality of the sample is further confirmed
that by a cortical fragility index of 1.26%.
[0411] Referring to FIG. 36, sample 2 is from a 65 year old woman.
The decay of the cortex is visible on the micrograph (top register)
and confirmed by image processor 16, which quantified CD1 and CD2
as 12% and 20%
[0412] respectively. The more decayed state of the sample is
further confirmed by the percentage of compact cortex (determined
by image processor 16 to be 57.1%), the percentage of
trabecularized cortex (PTC) (determined by image processor 16 to be
37.1%), and the index of cortical fragility (determined by image
processor 16 to be 6.8%).
[0413] In addition, DICOM files with scan data obtained from living
individuals (five from fracture cases and three from non-fracture
individuals) were obtained, to test system 10 with QCT data files
relating to living individuals. The relevant points in the density
profiles were identified and the relevant indices determined with
image processor 16. It was found that:
[0414] 1) Image processor 16 was able to identify relevant points
(i.e. j, k, l, m, n, o) and to compare image and direct SEM
measurements as discussed above.
[0415] Correlation between porosity measured directly and porosity
estimated by image processor 16 was excellent.
[0416] 2) The porosity (assessed by image processor 16) correlates
poorly with bone density, the currently used method to identify
individuals at risk for fracture. This suggests that the two
techniques overlap but do not capture the same things.
[0417] In the appendix (end of the text), are presented examples of
specimens studied. Micrographs and porosity scores as determined by
image processor 16 are shown as well as the BMD T-scores values in
these individuals. A BMD T-score of less than -1 is considered
osteopenic while T-score of -2.5 or less is osteoporosis.
[0418] Comparative values of porosity in individuals with fracture
of distal radius and those without fracture are shown in the
appendix also. The indices and parameters determined by image
processor 16 allow discrimination between individuals with and
without fracture-vulnerable bones.
[0419] In a further example, it was shown that the age-related
increase in porosity as measured using image processor 16 is
similar to that of activation frequency measured directly using
histomorphometry. FIG. 37A is a plot of age-related increase in
activation frequency (Ac.F) per year measured directly by
histomorphometry as a function of age (Compston, Porosity increase
from the rise in Ac.F within intracortical pores (haversian and
Volkmann canals, private communication). FIG. 37B is a plot of
age-related increase in porosity measured using image processor 16,
in 73 women. The data in these two figures demonstrates similar
trends.
[0420] FIGS. 38A, 38B and 38C are cross sectional views of a
radius. FIG. 38A is a greyscale image of the distal radius
reconstructed from high-resolution peripheral quantitative computed
tomography (HRpQCT) data. FIG. 38B is an image of the same sample,
reconstructed by image processor 16 after extraction of the
periosteal surface, demonstrating the ability of image processor 16
to segment an image from the background. The image was
reconstructed from BVE values; black (generally inner region) to
red (generally outer region) represent increasing BVE values. It
will be noted that, because a threshold is not used by image
processor 16, image processor 16 extracts the external contour of a
material (periosteal surface of the bone in this case) while
preserving the original attenuation of all the voxels within the
material. This allows image processor 16 to fully characterize the
material after it has been segmented from the background, including
assessing porosity and tissue mineralization level.
[0421] FIG. 38C is another view of the sample, showing the
segmentation of the same image into the four compartments of bone
described above. Red is the compact cortex 382, light green is the
outer transition zone 384, cyan is the inner transition zone 386,
and blue is the marrow space occupied by trabecular bone 388.
(Inner transition zone 386, where visible, lies between outer
transition zone 384 and trabecular bone 388. FIGS. 38D and 38E,
respectively, plot outer transition zone 384 and inner transition
zone 386 only in greyscale.) White is the background or surrounding
muscle tissue 390. FIGS. 38F and 38G are greyscale versions of,
respectively, FIGS. 38B and 38C.
[0422] In addition, an example of automated identification and
analysis of a foreign body within the bone is shown in FIGS. 39A,
39B and 40. A discussed above, image processor 16 can automatically
identify and analyze a structure within an image. While the main
application of system 10 is the identification and analysis of
bone, the use of image processor 16 is by no means limited to this
application.
[0423] FIG. 39A is an image as displayed by display 22 of image
processor 16. A bright object is evident within the image and,
although the bright object resembles a cross section of a bone, it
is actually the head of a metallic safety pin. This experiment was
designed to simulate the presence of a metallic foreign body within
a subject or other object, such as a bullet or a coin (perhaps
inadvertently swallowed).
[0424] The head of the safety pin was buried insight muscle tissue
and imaged with system 10; the bright object in the image of FIG.
39A is thus the cross section of the head of the safety pin. Image
processor 16 readily identified this structure, separated it from
the surrounding muscle tissue, and analyzed its structure,
permitting a diagnosis of a metallic body within the muscle
tissue.
[0425] ROI.sub.1--which is marked in FIG. 39A--is shown enlarged in
FIG. 39B. FIG. 40 (upper register) is the density profile curve
associated with ROI.sub.1. FIG. 40 (lower register) is the function
.lamda..sub.1 associated with the density profile curve shown in
the upper register of FIG. 40. As described above, .lamda..sub.1 is
the function that enables the bone to be distinguished from the
surrounding soft tissue. In this case .lamda..sub.1 enables the
separation of the object that was analyzed (viz. the head of the
safety pin) and the surrounding muscle tissue. The function
.lamda..sub.1 determines that a clear separation between the object
and the surrounding muscle is at point 99 (see the lower register
of FIG. 40). Point 99 is indicated in the density profile curve
(see the upper register of FIG. 40), and in the image of ROI.sub.1
(see FIG. 39B). It should be noted that j, in this case point 99,
is selected by .lamda..sub.1 in such a way that most of the
haziness (white blurry area adjacent to the object, open square)
adjacent to the structure is not included in the surrounding soft
tissue. This hazy area is due to partial volume effects, and image
processor 16--using .lamda..sub.1--can separate the object from the
surrounding muscle tissue with minimal artifacts due to partial
volume effect on muscle tissue. Furthermore, partial volume effects
do not affect the ability of image processor 16 to identify the
beginning of the object (point m, in this case at point 110) free
of artifacts. Image processor 16 readily identifies point n (at
point 115), enabling it to isolate the mass of the object free of
artifacts (PVE) and free decay (referred to as trabecularization in
the case of a bone).
[0426] Image processor 16 therefore identifies, and isolates the
object, and the surrounding muscle while minimizing the
interference from artifacts. The nature of the object can thus be
identified by comparing its density to that of the surrounding
muscle, whether by calculating its relative density or its absolute
density (such as in mgHA/cc or g/cc).
[0427] Image processor 16 determined the relative density of the
object viewed in FIG. 39A to be 93.5%. Hence, the object (in fact
the head of a safety pin) is much denser than muscle and is
therefore not bone (the relative density of which is typically
.about.66% and no greater than 75%, and .about.3 times greater than
muscle). In absolute terms, image processor determined the density
of the object to be 6036 mgHA/cc or 4.6 g/cc equivalents,
confirming that the object is of a density corresponding to a
metal. (In absolute terms, the density of purely mineralized bone
does not exceed .about.1200 mgHA/cc equivalent or .about.2 g/cc.)
The presence of an object of a much greater density than the bone
within a muscle tissue indicates the presence of a foreign
body.
[0428] The examples presented are not a complete list of the
potential used of image processor 16. It is envisaged that its
potential ability to automatically identify, separate and analyze
objects (or elements) of an unspecified image (i.e. where the
constituents are unknown), often without human intervention, will
have applications in other areas of medicine (such as in the
diagnosis of vessel calcifications, acute myocardial infraction,
fracture and stroke) but also outside the field of medicine. For
example, one such application is for the automated image analysis
of image processor 16 for the creation of artificial vision. In
these circumstances, image processor 16 would be integrated in a
camera to perform a contemporaneous analysis of images. After
analyzing the image, image processor 16 instantly transfers the
results to a "command centre" for decisions making. Such a "command
centre" could be a person or a software with specific instructions.
Artificial vision would have many applications within and outside
of medicine. Outside medicine, image processor 16 could be used,
for example, for controlling a robot.
[0429] The examples presented above suggest that image processor 16
is a useful tool for an automated identification, separation and
analyses of the structure of an object within an image, and in
particular the automated segmentation and analysis of bone. The
indices determined by image processor 16 allow the assessment of
cortical bone structure that is not the mere measurement of
cortical thickness, cortical area and cortical density.
[0430] Image processor 16 can retrieve an image from surrounding
background, while minimizing the confounding effects of PVE. For
example, the upper register of FIG. 41A is a cross sectional view
of a distal radius; the upper register of FIG. 41B is the same
cross section reconstructed by image processor 16 after separating
the bone (radius in this case) from the surrounding background. The
lower register of FIG. 41A is an enlarged view of an interface (or
junction) between the background and the bone, which is blurry
owing to partial volume effects (i.e., voxels that are partly bone,
partly background). The lower register of FIG. 41B, however, is the
same junction after processing of the image by image processor 16.
The delineation is much sharper: blurry voxels (i.e. those tainted
by partial volume effects) have been separated and left in the
background, and hence the margin of bone is more clearly
apparent.
[0431] This process of separation of a material within an image
from the background while minimizing the blurring due to partial
volume effects has many applications in medicine and outside the
medical field.
[0432] Image processor 16 can separate bone from background better
than existing software that is based on a threshold approach. FIG.
42A, upper register, is a cross sectional view of a tibia. The
external boundary is drawn with conventional threshold-based
software. It can be seen that a significant amount of bone has been
left out of this boundary and erroneously included in the
background. FIG. 42A, lower register, is a magnified view showing a
particular portion of bone outside the boundary and hence
erroneously included in the background.
[0433] FIG. 42B, upper register, is the same cross section of the
same bone separated from the background by image processor 16. It
can be seen that image processor 16 has clearly separated bone from
the surrounding soft tissue. The region erroneously segmented by
conventional threshold-based software is shown in the lower
register of this figure, in which it is apparent that image
processor 16 has correctly separated this region from the
background.
[0434] Modifications within the scope of the invention may be
readily effected by those skilled in the art. It is to be
understood, therefore, that this invention is not limited to the
particular embodiments described by way of example hereinabove.
[0435] In the claims that follow and in the preceding description
of the invention, except where the context requires otherwise owing
to express language or necessary implication, the word "comprise"
or variations such as "comprises" or "comprising" is used in an
inclusive sense, that is, to specify the presence of the stated
features but not to preclude the presence or addition of further
features in various embodiments of the invention.
[0436] Further, any reference herein to prior art is not intended
to imply that such prior art forms or formed a part of the common
general knowledge in any country.
* * * * *