# Patent

## Title
Method for assessing structural spatial frequencies using hybrid sampling with non-zero gradient for enhancement of selective sampling

## Abstract

The disclosed embodiments provide a method for acquiring MR data at resolutions down to tens of microns for application in in-vivo diagnosis and monitoring of pathology for which changes in fine tissue textures can be used as markers of disease onset and progression. Bone diseases, tumors, neurologic diseases, and diseases involving fibrotic growth and/or destruction are all target pathologies. Further the technique can be used in any biologic or physical system for which very high-resolution characterization of fine scale morphology is needed. The method provides rapid acquisition of selected values in k-space, with multiple successive acquisitions of individual k-values taken on a time scale on the order of microseconds, within a defined tissue volume, and subsequent combination of the multiple measurements in such a way as to maximize SNR. The reduced acquisition volume, and acquisition of only select values in k-space along selected directions, enables much higher in-vivo resolution than is obtainable with current MRI techniques.

## Background
REFERENCE TO RELATED APPLICATIONS
This application is a continuation in part of application Ser. No. 14/840,327 filed on Aug. 31, 2015 which relies on the priority of U.S. provisional application Ser. No. 62/044,321 filed on Sep. 1, 2014 entitled SELECTIVE SAMPLING MAGNETIC RESONANCE-BASED METHOD FOR ASSESSING STRUCTURAL SPATIAL FREQUENCIES, Ser. No. 62/064,206 filed on Oct. 15, 2014 having the same title and Ser. No. 62/107,465 filed on Jan. 25, 2015 entitled MICRO-TEXTURE CHARACTERIZATION BY MRI, the disclosures of which are incorporate herein by reference. The application additionally relies on the priority of provisional application Ser. No. 62/302,577 filed on Mar. 2, 2016 entitled METHOD FOR ASSESSING STRUCTURAL SPATIAL FREQUENCIES USING HYBRID SAMPLING WITH LOW OR INCREASED GRADIENT FOR ENHANCEMENT OF VERY LOW NOISE SELECTIVE SAMPLING WITH NO GRADIENT, the disclosure of which is incorporated herein by reference.


BACKGROUND
Field of the Invention
The herein claimed method relates to the field of diagnostic assessment of fine textures in biological systems for pathology assessment and disease diagnosis, and in material and structural evaluation in industry and in engineering research. More specifically, the invention employs a method for repeat measurement of k-values associated with the spatial organization of biologic tissue texture, with the MRI machine gradients turned off and for k-values in an narrowly associated neighborhood with a low gradient. This allows assessment of tissue texture on a time-scale on the order of a msec, whereby the problem of patient motion becomes negligible. The method enables in vivo assessment, towards diagnosis and monitoring, of disease and therapy-induced textural changes in tissue. Representative targets of the technique are: 1) for assessment of changes to trabecular architecture caused by bone disease, allowing assessment of bone health and fracture risk, 2) evaluation of fibrotic development in soft tissue diseases such as, for example, liver, lung, and heart disease, and 3) changes to fine structures in neurologic diseases, such as the various forms of dementia, or in cases of brain injury and downstream neuro-pathology as in, for example, Traumatic Brain Injury (TBI) and Chronic Traumatic Encephalopathy (CTE), or for characterization and monitoring of abnormal neurologic conditions such as autism and schizophrenia. Other pathology applications include assessment of vascular changes such as in the vessel network surrounding tumors or associated with development of CVD (Cerebrovascular Disease), and of changes in mammary ducting in response to tumor growth. The invention also has applications in assessment of fine structures for a range of industrial purposes such as measurement of material properties in manufacturing or in geology to characterize various types of rock, as well as other uses for which measurement of fine structures/textures is needed.
Description of the Related Art
Though fine textural changes in tissue have long been recognized as the earliest markers in a wide range of diseases, robust clinical assessment of fine texture remains elusive, the main difficulty arising from blurring caused by subject motion over the time required for data acquisition.
Early and accurate diagnosis is key to successful disease management. Though clinical imaging provides much information on pathology, many of the tissue changes that occur as a result of disease onset and progression, or as a result of therapy, are on an extremely fine scale, often down to tens of microns. Changes in fine tissue texture have been recognized for many years by diagnosticians, including radiologists and pathologists as the earliest harbinger of a large range of diseases, but in vivo assessment and measurement of fine texture has remained outside the capabilities of current imaging technologies. For instance, differential diagnosis of obstructive lung disease relies on a textural presentation in the lung parenchyma, but the robustness of the Computed Tomography (CT) measure of early stage disease is limited. Trabecular bone microarchitecture, the determinant of fracture risk in aging bone, has also remained elusive due to image blurring from patient motion during Magnetic Resonance (MR) imaging scans. Post processing analysis of MR-images is sometimes used to try to differentiate surface textures in structures such as tumors and white matter. (DRABYCZ, S., et al.; “Image texture characterization using the discrete orthogonal S-transform”; Journal of Digital Imaging, Vol. 22, No 6, 2009. KHIDER, M., et al.; “Classification of trabecular bone texture from MRI and CT scan images by multi-resolution analysis”; 29th Annual International Conference of the IEEE Engineering in Medicine and Biology Society, EMBS 2007.) But post processing analysis is limited in effect as it doesn't deal with the underlying problem that prevents high resolution acquisition of textural information, i.e. subject motion. (MACLAREN, J. et al.; “Measurement and correction of microscopic head motion during magnetic resonance imaging of the brain”, PLOS/ONE, Nov. 7, 2012. MACLARAN, J. et al.; “Prospective motion correction in brain imaging: a review; Magnetic Resonance in Medicine, Vol. 69, 2013.)
The main sources of motion affecting MR imaging are cardiac pulsatile motion, respiratory-induced motion and twitching. The first two are quasi-cyclic, the usual approach to which is gating at the slowest phase of motion. However, even with gating, there is sufficient variation between acquisitions to cause loss of spatial phase coherence at the high k-values of interest for texture measurements. This problem is exacerbated by the fact that motion may not be perfectly cyclic, and often originates from combined sources. Twitching is rapid, inducing random displacements, and hence it is not possible to maintain coherence at the high k-values of interest when measuring texture.
While Positron Emission Tomography (PET) provides valuable diagnostic information, it is not capable of resolution below about 5 mm and relies on the use of radioactive tracers for imaging as well as x-ray beams for positioning, raising dose concerns, especially if repeat scanning is needed. (BERRINGTON DE GONZALEZ, A. et al.; “Projected cancer risks from Computed Tomographic scans performed in the United States in 2007”; JAMA Internal Medicine, Vol. 169, No. 22, December 2009.) Further, PET imaging is extremely costly, requiring a nearby cyclotron. CT resolution down to 0.7 mm is possible in theory, though this is obtained at high radiation dose and is subject to reduction by patient motion over the few minute scan time. The non-negligible risk from the associated radiation dose makes CT problematic for longitudinal imaging and limits available resolution. Along with serious dose concerns, digital x-ray resolution is limited because the 2-dimensional image obtained is a composite of the absorption through the entire thickness of tissue presented to the beam. Current clinical diagnostics for the diseases that are the target of the method claimed herein are fraught with difficulties in obtaining sufficient in vivo resolution, or accuracy. In some cases, no definitive diagnostic exists currently. In other pathologies, particularly in breast and liver, diagnosis is dependent on biopsy, with its non-negligible risk of morbidity and even mortality, and which is prone to high read and sampling errors. (WELLER, C; “Cancer detection with MRI as effective as PET-CT scan but with zero radiation risks”; Medical Daily, Feb. 18, 2014.)
Bone health is compromised by ageing, by bone cancer, as a side effect of cancer treatments, diabetes, rheumatoid arthritis, and as a result of inadequate nutrition, among other causes. Bone disease affects over ten million people annually in the US alone, adversely affecting their quality of life and reducing life expectancy. For assessment of bone health, the current diagnostic standard is Bone Mineral Density (BMD), as measured by the Dual Energy X-ray Absorptiometry (DEXA) projection technique. This modality yields an areal bone density integrating the attenuation from both cortical and trabecular bone, similar to the imaging mechanism of standard x-ray, but provides only limited information on trabecular architecture within the bone, which is the marker linked most closely to bone strength. (KANIS, J. AND GLUER, C.; “An update on the diagnosis and assessment of osteoporosis with densitometry”; Osteoporosis International, Vol. 11, issue 3, 2000. LEGRAND, E. et al.; “Trabecular bone microarchitecture, bone mineral density, and vertebral fractures in male osteoporosis”; JBMR, Vol. 15, issue 1, 2000.) BMD correlates only loosely with fracture risk. A post-processing technique, TBS (Trabecular Bone Score) attempts to correlate the pixel gray-level variations in the DEXA image, to yield information on bone microarchitecture. A comparison study determined that BMD at hip remains a better predictor of fracture. But, though TBS does not yield a detailed assessment of trabecular architecture. (BOUSSON, V., et al.; “Trabecular Bone Score (TBS): available knowledge, clinical relevance, and future prospects”; Osteoporosis International, Vol. 23, 2012. DEL RIO, et al.; “Is bone microarchitecture status of the spine assessed by TBS related to femoral neck fracture? A Spanish case-control study”: Osteoporosis International, Vol. 24, 2013.) TBS is a relatively new technique and is still being evaluated.
Measurement of bone microarchitecture, specifically trabecular spacing and trabecular element thickness, requires resolution on the order of tenths of a millimeter. MRI, ultrasound imaging, CT, and microCT have all been applied to this problem. In MRI, though high contrast between bone and marrow is readily obtained, resolution is limited by patient motion over the long time needed to acquire an image with sufficient resolution to characterize the trabecular network. The finer the texture size of this network, the greater the blurring from motion. An attempt to mitigate the effects of patient motion by looking only at the skeletal extremities, removed from the source of cardiac and respiratory motion sources, has been tried using both MRI and microCT. However, the correlation between bone microarchitecture in the extremities and that in central sites in not known. Further, a large data matrix, hence long acquisition time, is still required to obtain sufficient image information to determine trabecular spacing and element thickness. This long acquisition time results in varying levels of motion-induced blurring, depending on patient compliance—twitching is still a serious problem even when measuring extremities. A proposed MR-based technique, fineSA (JAMES, T., CHASE, D.; “Magnetic field gradient structure characteristic assessment using one dimensional (1D) spatial-frequency distribution analysis”; U.S. Pat. No. 7,932,720 B2; Apr. 26, 2011.), attempts to circumvent the problem of patient motion by acquiring a much smaller data matrix of successive, finely-sampled, one-dimensional, frequency-encoded acquisitions which are subsequently combined to reduce noise. Imaging in this case is reduced to one dimension, reducing the size of the data matrix acquired and, hence, the acquisition time. However, as the gradient encoded echoes, are very low Signal to Noise (SNR), noise averaging is required. Though some resolution advantage is gained by this method relative to 2 and 3-d imaging, the need to acquire many repeat spatially-encoded echoes over several response times (TRs) for signal averaging results in an acquisition time on the order of minutes—too long to provide motion immunity. Thus, resolution improvement obtainable by the technique is limited.
What is needed is an accurate, robust, non-invasive, in vivo measure of trabecular spacing and trabecular element thickness capable of assessing bones in the central skeleton, as these are the key markers for assessing bone health and predicting fracture risk. Until now, no clinical technique has been able to provide this capability.
Fibrotic diseases occur in response to a wide range of biological insults and injury in internal organs, the development of collagen fibers being the body's healing response. The more advanced a fibrotic disease, the higher the density of fibers in the diseased organ. Fibrotic pathology occurs in a large number of diseases, from lung and liver fibrosis, to cardiac and cystic fibrosis, pancreatic fibrosis, muscular dystrophy, bladder and heart diseases, and myelofibrosis, in which fibrotic structures replace bone marrow. Fibrotic development is attendant in several cancers, such as breast cancer. A different pathology development is seen in prostate cancer, where the disease destroys healthy organized fibrous tissue. In all cases, textural spacings highlighted in the tissue change in response to disease progression, as collagen fibers form along underlying tissue structures. In liver disease, the textural wavelength changes as the healthy tissue texture in the liver is replaced by a longer wavelength texture originating from the collagen “decoration” of the lobular structure in the organ. In other organs/diseases, textural change reflects the upset in healthy tissue with development of texture indicative of fibrotic intervention.
To span the range of disease progression in most fibrotic pathologies, evaluation of textural changes from fibrotic development requires resolution on the scale of tenths of a mm One of the most prevalent of such pathologies, liver disease, is representative of the difficulty of assessing fibrotic structure. Currently, the gold standard for pathology assessment is tissue biopsy—a highly invasive and often painful procedure with a non-negligible morbidity—and mortality—risk (patients need to stay at the hospital for post-biopsy observation for hours to overnight), and one that is prone to sampling errors and large reading variation. (REGEV, A.; “Sampling error and intraobserver variation in liver biopsy in patients with chronic HCV infection”; American Journal of Gastroenterology; 97, 2002. BEDOSSA, P. et al.; “Sampling variability of liver fibrosis in chronic hepatitis C”; Hepatology, Vol. 38, issue 6, 2004. VAN THIEL, D. et al.; “Liver biopsy: Its safety and complications as seen at a liver transplant center”; Transplantation, May 1993.) Ultrasound, another modality often used to assess tissue damage in liver disease, is only able to provide adequate assessment in the later stages of the disease—it is used to diagnose cirrhosis. Magnetic Resonance-based Elastography (MRE), which has been under development for some time for use in assessment of liver disease, is not capable of early-stage assessment—the read errors are too large prior to significant fibrotic invasion (advanced disease). Further, this technique requires expensive additional hardware, the presence of a skilled technician, and takes as much as 20 minutes total set up and scanning time, making it a very costly procedure. The ability to image fibrotic texture directly by MR imaging is compromised both by patient motion over the time necessary to acquire data and by lack of contrast between the fibers and the surrounding tissue. Even acquisition during a single breath hold is severely compromised by cardiac pulsatile motion and noncompliance to breath hold, which results in significant motion at many organs, such as liver and lungs. And SNR is low enough that motion correction by combining reregistered MR-intensity profiles obtained from successive echoes is extremely problematic. Similarly, assessment of the amount of cardiac fibrosis in early stage disease using MRI is seriously hampered by cardiac pulsation over the time of the measurement. As motion is, unlike Gaussian noise, a non-linear effect, it can't be averaged out—there must be sufficient signal level to allow reregistration before averaging for electronic noise-reduction. A more sensitive (higher SNR), non-invasive technique, capable of assessing textural changes throughout the range of fibrotic development, from onset to advanced pathology, is needed to enable diagnosis and monitoring of therapy response.
Onset and progression of a large number of neurologic diseases are associated with changes in repetitive fine neuronal and vascular structures/textures. However, ability to assess such changes in the brain is only available post mortem. Currently, definitive diagnosis of Alzheimer's Disease (AD) is by post mortem histology of brain tissue. AD and other forms of dementia such as Dementia with Lewy Bodies, motor diseases such as Amyotrophic Lateral Sclerosis (ALS), Parkinson's disease, conditions precipitated by Traumatic Brain Injury (TBI) such as Chronic Traumatic Encephalopathy (CTE), as well as those caused by other pathologies or trauma, or conditions that involve damage to brain structures such as Multiple Sclerosis (MS), Cerebrovascular Disease (CVD), and other neurologic diseases, are often only diagnosable in advanced stages by behavioral and memory changes, precluding the ability for early stage intervention. Further, conditions such as epilepsy and autism have been associated with abnormal variations in fine neuronal structures, which, if clinically diagnosable, would allow targeted selection for testing therapy response.
Various in vivo diagnostic techniques are available for AD and other dementias, but none of them are definitive. These techniques range from written diagnostic tests, which are prone to large assessment errors, to PET imaging to assess amyloid plaque density or glucose metabolism (FDG PET). As discussed previously, PET imaging is extremely expensive, cannot provide high resolution, and relies on use of radioisotopes and positioning x-ray beams, complicating approval for longitudinal use due to dose concerns. Further, neither amyloid imaging nor FDG PET has been shown to provide a definitive indication of AD. (MOGHBEL, M. et al. “Amyloid Beta imaging with PET in Alzheimer's disease: is it feasible with current radiotracers and technologies?”; Eur. J. Nucl. Med. Mol. Imaging.)
Use of CSF biomarkers for dementia diagnosis is painful and highly invasive and cannot differentiate signal levels by anatomic position in the brain, as is possible with imaging biomarkers. As various forms of dementia are found to have different spatial/temporal progression through the brain, this is a serious drawback to use of liquid biopsy. Another disease associated with various forms of dementia is CVD (Cerebrovascular Disease), which induces cognitive impairment as a result of reduced blood flow through blocked vessels leading to brain tissue. Something capable of high-resolution assessment of pathology-induced changes in micro-vessels is needed here.
Tissue shrinkage due to atrophy in many forms of dementia including AD is measurable with careful registration of longitudinally-acquired data in MRI, but the disease is advanced by the time this shrinkage is measureable. Early stages of disease are indicated in post mortem histology by degradation in the columnar ordering of cortical neurons, the normal spacing for these columns being on the order of 100 microns in most cortical regions. (CHANCE, S. et al.; “Microanatomical correlates of cognitive ability and decline: normal ageing, MCI, and Alzheimer's disease”; Cerebral Cortex, August 2011. E. DI ROSA et al.; “Axon bundle spacing in the anterior cingulate cortex of the human brain”; Journal of Clinical Neuroscience, 15, 2008.) This textural size, and the fact that the cortex is extremely thin, makes speed of acquisition paramount, as even tiny patient motion will make data collection impossible. Assessment of textural changes on the order of tens of microns microns is extremely problematic in vivo, but would, if possible, enable targeting a range of fine textural changes in neuronal disease diagnosis and monitoring, and would play an important role in therapy development.
Another possible neurologic application for the claimed method is to, in vivo, determine the boundaries of the various control regions of the cerebral cortex or the different Brodmann's areas of which these are comprised. Such ability would greatly aid data interpretation in brain function studies, such as those performed using, for example, FMRI (Functional Magnetic Resonance Imaging).
The three classes of diseases listed above, bone disease, fibrotic diseases, and neurologic diseases are not an all-inclusive list. Other disease states in which pathology-induced changes of fine structures occur, for instance angiogenic growth of vasculature surrounding a tumor, or fibrotic development and changes in vasculature and mammary gland ducting in response to breast tumor development, also are pathologies wherein the ability to resolve fine tissue textures would enable early detection of disease, and monitoring of response to therapy.
The ability to measure changes in fine textures would be of great value for disease diagnosis. Non-invasive techniques that do not rely on use of ionizing radiation or radioactive tracers allow the most leeway for early diagnosis and repeat measurement to monitor disease progression and response to therapy. Magnetic Resonance Imaging (MRI), which provides tunable tissue contrast, is just such a non-invasive technique, with no radiation dose concerns. However, in order to circumvent the problem of signal degradation due to patient motion, data must be taken on a time scale not previously possible.
SUMMARY OF THE INVENTION
A method for selective sampling to assess texture using magnetic resonance (MR) is accomplished by exciting a volume of interest. All gradients are then turned off and multiple samples of an RF signal encoded at a specific k-value are recorded for a very low SNR region. Hybrid addition of measurement at a first low SNR region of k-values in a narrow neighborhood using a small non-zero gradient is then employed to enhance the data. Further addition of measurements at a high SNR region using a higher non-zero gradient may then additionally employed.



BRIEF DESCRIPTION OF THE DRAWINGS
The features and advantages of embodiments disclosed herein will be better understood by reference to the following detailed description when considered in connection with the accompanying drawings wherein:
FIG. 1 is a simulation showing the number of data samples required for averaging to achieve an output SNR≧20 dB as a function of input SNR;
FIG. 2 is a simulation showing the number of data samples needed for averaging to achieve a SNR≧20 db as a function of location in k-space;
FIG. 3 is an example timing diagram of a pulse sequence for the claimed method showing the timing of a single TR;
FIG. 4 is a close-up of the example timing diagram of FIG. 3;
FIG. 5 is an example of a timing diagram for the claimed method, designed to acquire multiple measures of a select set of k-values, with a different number of samples acquired at each k-value to counteract the decrease in energy density at increasing k-value;
FIG. 6 is a simulation showing that the ability provided by the claimed method to acquire many repeats of a targeted k-value within a single TR enables robust signal averaging to boost SNR;
FIG. 7 is a simulation showing the results of attempting to acquire 90 samples for averaging using the conventional frequency-encoded echo approach, wherein acquisition of only a small number of repeats of a particular k-value are possible in each TR due to the long record time for each echo;
FIG. 8 is an example timing diagram for the claimed method designed to provide data acquisition over multiple refocused echoes within a single TR; and,
FIGS. 9 and 10 are a depiction of two possible shapes for the acquisition volume of interest (VOI);
FIG. 11 is an example timing diagram of a pulse sequence for the claimed hybrid method showing the timing of a single TR;
FIG. 12 is a detailed view of the hybrid elements of the method at an expanded scale;
FIG. 13 is a further detailed view of the very-low SNR acquisition mode portion of FIG. 12;
FIG. 14 is a further detailed view of the low SNR acquisition portion of FIG. 12;
FIG. 15 is a further detailed view of the high SNR acquisition portion of FIG. 12;
FIG. 16 is an example timing diagram of a pulse sequence for the claimed hybrid method showing data acquisition in a single echo;
FIG. 17 is a further detailed view of the very-low SNR, low SNR and high SNR acquisition portions of FIG. 16;
FIG. 18 is a an example timing diagram of a pulse sequence for a low SNR acquisition;
FIG. 19 is a further detailed view of the low SNR acquisition mode of FIG. 18;
FIG. 20 is an example timing diagram of a pulse sequence for high SNR acquisition; and,
FIG. 21 is a further detailed view of the high SNR acquisition mode of FIG. 20.



DETAILED DESCRIPTION OF THE INVENTION
The following definition of terms as used herein is provided:
180° inversion pulse RF pulse that inverts the spins in a tissue region to allow refocusing of the MR signal.
180° pulse An RF pulse that tips the net magnetic field vector antiparallel to B0 
90° pulse An RF pulse that tips the net magnetic field vector into the transverse plane relative to B0 
3 T 3 Tesla
A/D Analog to digital converter
AD Alzheimer's Disease
ADC Average diffusion coefficient measured in Diffusion Weighted Imaging
Adiabatic pulse excitation Adiabatic pulses are a class of amplitude- and frequency-modulated RF-pulses that are relatively insensitive to B1 inhomogeneity and frequency offset effects.
AWGN Additive White Gaussian Noise Additive white Gaussian noise (AWGN) is a basic noise model used in Information theory to mimic the effect of many random processes that occur in nature.
Biopsy A biopsy is a sample of tissue extracted from the body in order to examine it more closely.
C/N Contrast to Noise, a measure of image quality based on signal differences between structural elements rather than on overall signal level
CAWGN Complex-valued, additive white Gaussian noise
CBF Cerebral Blood Flow
Chemical shift Small variations in MR resonant frequency due to the different molecular environments of the nuclei contributing to an MR signal.
CJD Creutzfeld-Jakob Disease
Crusher gradients Gradients applied on either side of a 180° RF refocussing slice selection pulse to reduce spurious signals generated by imperfections in the pulse.
CSF Cerebrospinal fluid
DEXA Dual Energy X-ray Absorptiometry is a means of measuring bone mineral density using two different energy x-ray beams.
DSC Dynamic susceptibility contrast
DTI Diffusion Tensor Imaging
DWI Diffusion Weighted Imaging
Echo The RF pulse sequence where a 90° excitation pulse is followed by a 180° refocusing pulse to eliminate field inhomogeneity and chemical shift effects at the echo.
Frequency encodes Frequency-encoding of spatial position in MRI is accomplished through the use of supplemental magnetic fields induced by the machine gradient coils
Gaussian noise Gaussian noise is statistical noise having a probability density function (PDF) equal to that of the normal distribution, which is also known as the Gaussian distribution.
Gradient pulse a pulsing of the machine magnetic field gradients to alter the k-value encode
Gradient set the set of coils around the bore of an MR scanner used primarily to spatially encode signal or to set a particular phase wrap in a selected direction
GRE Gradient Recalled Echo
Interleaved acquisition Signal acquisition from a multiplicity of VOIs, successively excited within a single TR
Isochromat A microscopic group of spins that resonate at the same frequency.
k-space The 2D or 3D Fourier transform of the MR image.
k-value coefficient The coefficient in a Fourier series or transform reflecting the relative weight of each specific k-value in the series.
k-space The 2D or 3D Fourier transform of the MR image.
k-value One of the points in k-space reflecting the spacing of structural elements in a texture field.
k-value selection pulse The gradient pulse used to select a specific k-value encode along the sampled direction
Library of k-space values the net collection of k-space coefficients acquired in a particular region of tissue for tissue characterization
Machine gradients the magnetic field gradients imposable through use of the set of gradient coils in an MR scanner
MRE Magnetic Resonance Elastography—an imaging technique that measures the stiffness of soft tissues using acoustic shear waves and imaging their propagation using MRI.
MRI Magnetic Resonance Imaging
MS Multiple Sclerosis
Noise floor In signal theory, the noise floor is the measure of the signal created from the sum of all the noise sources and unwanted signals within a measurement system
PET Positron Emission Tomography is a functional imaging technique that produces a three-dimensional image of functional processes in the body using a positron-emitting radiotracer.
Phase coherence (spatial) When referring to multiple measurements within a common VOI of a or multiple k-values indicates that the sample has the same position relative to the measurement frame of reference
Phase encode A phase encode is used to impart a specific phase angle to a transverse magnetization vector. The specific phase angle depends on the location of the transverse magnetization vector within the phase encoding gradient, the magnitude of the gradient, and the duration of the gradient application.

Phase wrap The helical precession of the phase of the transverse magnetization along a phase encoded sample

Pitch with reference to the pitch of a screw, the tightness of the phase wrap along the direction of k-value encode

Profile A one dimensional plot of signal intensity

RF Radio Frequency electromagnetic signal

Semi-crystalline texture a texture exhibiting regular spacing along one or more directions

slice (slab) Used interchangeably to indicate a non-zero thickness planar section of the

Slice-selective refocusing Refocussing of spins through combination of a slice selective gradient and an RF pulse such that the bandwidth of the RF pulse selects a thickness along the direction of the gradient, and the RF pulse tips the net magnetization vector away from its equilibrium position Only those spins processing at the same frequency as the RF pulse will be affected.

SE Spin Echo

SNR Signal to Noise Ratio

Spoiler gradients see crusher gradients

T2 Defined as a time constant for the decay of transverse magnetization arising from natural interactions at the atomic or molecular levels.

T2* In any real NMR experiment, the transverse magnetization decays much faster than would be predicted by natural atomic and molecular mechanisms; this rate is denoted T2* (“T2-star”). T2* can be considered an “observed” or “effective” T2, whereas the first T2 can be considered the “natural” or “true” T2 of the tissue being imaged. T2* is always less than or equal to T2.

TBS Trabecular Bone Score is a technique that looks for texture patterns in the DEXA signal for correlation with bone microarchitecture for assessing bone health

TE Spin Echo sequences have two parameters: Echo Time (TE) is the time between the 90° RF pulse and MR signal sampling, corresponding to maximum of echo. The 180° RF pulse is applied at time TE/2. Repetition Time is the time between 2 excitations pulses (time between two 90° RF pulses).

Textural frequency the number of texture wavelength repeats per unit length in a texture

Texture wavelength the characteristic spacing between structural elements in a texture

TR Spin Echo sequences have two parameters: Echo Time (TE) is the time between the 90° RF pulse and MR signal sampling, corresponding to maximum of echo. The 180° RF pulse is applied at time TE/2. Repetition Time is the time between 2 excitations pulses (time between two 90° RF pulses).

Vector combination gradient A magnetic gradient resulting from any vector combination of the gradient coil set

VOI Volume of Interest

Windowing function In signal processing, a window function (also known as an apodization function or tapering function) is a mathematical function that is zero-valued outside of some chosen interval

x-ray diffraction X-ray diffraction is a tool used for identifying the atomic and molecular structure of a crystal

The embodiments disclosed herein provide an MR-based technique that enables in vivo, non-invasive, high-resolution measurement and assessment of fine biologic textures, enabling monitoring of texture formation and/or change in response to disease onset and progression in a range of pathologies. This same method can be applied to fine-texture characterization in other biologic and physical systems. It enables MR-based resolution of fine textures to a size scale previously unattainable in in vivo imaging. The method, while described herein with respect to biological systems for examination of tissue, is equally applicable for assessment of fine structures in a range of industrial purposes such as measurement of material properties in manufacturing or in geology to characterize various types of rock, as well as other uses for which measurement of fine structures/textures is needed.
The method claimed herein achieves this significant improvement in in vivo resolution of fine texture by acquiring the requisite data fast enough that the effect of subject motion, the factor that limits MRI resolution, becomes negligible. This fast acquisition is achieved by acquiring data incrementally—at a single location, orientation and at one, or a select set, of k-values at a time—within one TR. After applying an encoding gradient to select the k-value of interest, data is acquired with the gradient switched off, allowing multiple acquisition repeats of the encoded k-value for subsequent averaging to reduce electronic noise, thus enabling robust measure of individual k-values before motion blurring can occur. To build up measurements on a larger set of selected k-values present within the tissue, or towards development of a continuous spectrum of textural spacings within the tissue, the acquisition TR can be repeated as many times as necessary, changing the encode as needed to span the desired extent of real and of k-space required. The set of one or more k-values output from each TR are now high SNR due to the ability to average repeats without motion effects, and since the measure of interest is textural spacing, and not development of an image, the lack of phase coherence between TRs is of no concern.
In its simplest form the method claimed herein consists of acquiring MR signal from within an inner volume to encompass a specific tissue region of interest, such as a lesion, an organ, a location in an organ, a specific region of bone, or a number of regions in a diseased organ for sampling. This inner volume may be excited by one of a number of methods, including but not limited to: intersecting slice-selective refocusing, selective excitation using phased-array transmit in combination with appropriate gradients, adiabatic pulse excitation to scramble signal from the tissue outside the region of interest, outer volume suppression sequences, and other methods of selectively exciting spins in an internal volume including physically isolating the tissue of interest, to name a few,
After definition of a volume of interest (VOI), the gradient is turned off, and multiple samples of signal centered at a specific k-value, the spread of which is defined by receiver BW and sampling length, are acquired. This measurement is repeated only in specified directions within the VOI rather than trying to map all of k-space, as is required to generate an image. One or more samples of a particular k-value are acquired within an acquisition block during a single TR and the k-value subsequently incremented or decremented, allowing further multiple samples of other k-values as desired during the same TR. This method allows multiple sampling of each k-value of interest over a time period of milliseconds, providing immunity to subject motion. The process can then be repeated in further TRs, the requirement on motion between k-value acquisitions being only that the VOI remain within the tissue region of interest. Build up of a magnitude spectrum of spatial frequencies may be accomplished without the need to acquire it in a spatially coherent manner. Because the quantities of interest are the relative intensities of the various k-values (textural spacings) present in the sample volume, as long as the acquisition volume remains within a representative sample of tissue, any motion between the blocks does not compromise the measurement. In the case where motion of sufficiently large magnitude that the internally excited volume could move into other tissue volumes over the course of building up a spectrum of k-values contained in the tissue, use of fairly robust, real-time piloting and acquisition algorithms can be used for gross repositioning of the internal selectively excited volume and for rejecting data sets that have failed to stay in the proper tissue.
Repositioning the VOI to allow sampling of texture at multiple positions within or across an organ or anatomy allows determination of the variation in pathology through the organ. The data acquired can, with reference to positioning images, be mapped spatially. Either the VOI can be moved in successive TRs or interleaved acquisition done within a single TR by exciting additional volumes during the time that the signal is recovering in advance of the next TR. The requirement is that successive VOIs be excited in new tissue, that does not overlap the previous slice selects. Spatial variation of pathology can be determined by this method. This can also be used to monitor temporal progression of a pathology through an organ if the measure is repeated longitudinally.
Tailoring the pulse sequence to pre-wind phase in the sample volume positions the highest k-values of interest at the echo peak where the signal is strongest, providing best SNR measurement.
Sampling of k-values along multiple directions at varying angles and along varying paths, either rectilinear or curved, within the volume under study can yield important information on texture, especially textures with semi-ordered structure in specific directions, such as neuronal minicolumns. Measurement of the k-values associated with columnar spacing is extremely sensitive to alignment of the sampling path, as slight variations in sampling direction on either side of perpendicular show a rapid drop off in signal magnitude for that k-value. Rocking the acquisition path on either side of the signal maximum can yield a measure of pathology-induced randomness which is indicated by the width of the peak.
With the gradient switched off for data acquisition, tuning the bandwidth to particular chemical species can enhance structural information when the chemical composition of the structure under study is known.
The method claimed herein can be used in conjunction with time-dependent contrast schemes that target blood flow. Some of these contrast techniques are Blood Oxygenation Level Dependent (BOLD) imaging, Arterial Spin Labeling (ASL) imaging, and Dynamic Susceptibility Contrast (DSC) imaging. As these methods use various techniques to highlight vasculature, changes in the texture of the vasculature associated with many pathologies, including CVD (cerebrovascular disease) and tumor growth can be measured.
The method claimed herein can also be used in conjunction with other MR-based measurement techniques, including DWI and DTI, to provide front end information toward parameter selection for the diffusion techniques as well as correlation with their measurements of tissue health.
The rapid repeat measurement of a single k-value, with the total time to acquire a block being on the order of a msec, reduces patient and machine motion-induced blurring to a negligible level, enabling robust assessment of fine textures previously not accessible in vivo. (For comparison, standard MR image acquisition times are much longer in duration over which patients are asked to remain completely stationary.) The SNR of each k-value measured is significantly improved through combination of the individual samples at each k-value within a block; this averaging can now can be done without concern for subject motion, which is eliminated due to the rapid sequential acquisition of the individual samples in the block.
This significant improvement in SNR is made possible because the method claimed herein focuses on acquiring only the k-values of interest for determination of fine texture pathology signatures rather than on acquisition of the large number of entire spatially-encoded echoes required for image formation. The significantly reduced data matrix enables the increased number of repeats at the targeted k-values, and hence significant improvement in SNR.
Energy density within a range of textural spacings is proportional to textural wavelength, or inversely to k-value—i.e. the higher the k-value, the lower the associated signal intensity. The fast acquisition enabled through use of the method claimed herein, enables tailoring the number of acquisition repetitions at a particular k-value to acquire k-values for which there is low signal first, before T2 and T2* effects have degraded signal amplitude. In this way, the SNR of each repeat to be averaged for noise cancellation (or spatial-phase-corrected before combining it with the measurements of k-value from subsequent TRs) will be above this threshold. It does not matter that there is motion between acquisition cycles at different k-values as long as each acquisition lies within the tissue volume of interest (VOI). As the claimed method targets only assessment of pathology-induced changes in tissue texture, there is no requirement for phase coherence over an entire cycle of data acquisition, as is required in imaging.
Several benefits result from acquiring data after the gradient is switched off for single-k-value sampling in a reduced volume (the VOI). By proper pulse sequencing, the echo record window can be designed such that recording begins with the highest k-values of interest, as signal level is highest at the echo peak. This enables recording of fine structures currently unachievable with in vivo MR imaging.
Additionally, T2* is longer with the gradient off, so SNR is improved by the longer acquisition times possible This allows acquisition of an increased number of samples, N.
Coil combination is also simplified by having higher SNR for each k-value, hence providing a significant improvement in overall SNR. This is especially beneficial as the trend in MRI is towards coil arrays composed of many small element coils. As the acquisition volumes targeted in the method claimed herein are small, correction for phase across the sample volume is not needed. Only one phase and gain value for each coil is needed for combining the multiple element channels. These can be combined using the Maximal Ratio Combining (MRC) method, which weights the coil with the highest SNR most heavily, or other multi-signal combination methods. (Phase and gain for the elements of a given coil array can be determined once from a phantom and applied to patient data.)
Signal acquisition and data sampling in a standard MRI scan is done by acquiring complex-valued samples of multiple echoes, while applying a gradient sequence concurrently, as well as in sequence with the echoes. Imaging relies on frequency encode for one of the dimensions because this allows a line in k-space to be acquired with each phase encode rather than a single point. For 3-dimensional imaging, two dimensions in k-space normally rely on phase encode to generate the targeted filling of k-space, with the third dimension frequency-encoded. Phase encode acquisition in imaging usually entails acquisition of on the order of 256 k-values in each of the phase-encode directions, hence is is a relatively slow process. Clinical MRI scans take on the order of 10-15 minutes to generate an image. The aim in image construction is to acquire sufficient k-space coverage to fill out all the coefficients in the 2 or 3-dimensional Fourier series, which is why in standard MR resolution is limited by subject motion.
The method claimed herein is in direct contrast to standard MR data acquisition, with its focus on image generation. Image formation is plagued by blurring resulting from subject motion over the long time necessary to acquire the large data matrix required. Since the target of the method claimed herein is texture rather than image, the only requirement on subject motion is that the sampled volume remain within a region of similar tissue properties over the course of acquiring data. This is a much less stringent and easy to achieve target than the requirement of structural phase coherence, as the scale of the allowable motion is then large enough, and of a temporal order, to be easily correctable by real-time motion assessment and correction techniques. The speed of acquisition for the method claimed herein is such that, in most cases, real-time motion correction may not be necessary at all. While other methods have focused on post-processing of images to try to extract textural measures, the method claimed herein eliminates the need for image generation, focusing instead on directly measuring texture, hence enabling a more sensitive and robust measure.
Frequently, k-space sampling is considered synonymous with sampling of an echo in the presence of a gradient set. In the method claimed herein, the approach to k-space filling is to acquire only the set of k-values needed for texture evaluation in the targeted pathology, with data acquired after the gradient is switched off. This method enables such rapid acquisition of single-k-value repeats for averaging for noise reduction that subject motion does not degrade the data.
Along with the huge improvement in SNR that arises from sampling k-values individually, with many repeats of a select set of k-values acquired in a single TR, acquisition after the gradient is switched off allows further significant improvement in SNR and hence, increase in measurement robustness. This is explained in the following discussion.
MR echo sampling provides specific samples vs. time of a time-dependent echo. The echo is comprised by the gradients applied concurrently (for the frequency-encode axis) and prior to (for a phase encode axis), but also contains the isochromats associated with the different chemical species of the sample, as well as the envelope (T2 & T2*) associated with spin-spin interactions.
Conventional frequency-encoded spin acquisitions impose a time-varying gradient upon the sample, which effectively travels in k-space along a pre-defined path. For rectilinear sampling, the path is along a straight line.
Frequency encodes generate only one measurement at a given k-value—at a given point in time, the acquired sample of the echo represents the one value which corresponds to the Fourier coefficient at a specific k-space location. The next echo sample represents the value at a different k-space location, the next k-value dependent on the slope of the gradient applied concurrently. As long as there is sufficient signal at the corresponding k-value, this approach works well. However, in cases where the signal of interest is near or even below the noise floor, usually additional samples and subsequent post-processing will be required.
One way to reduce the noise floor in a frequency-encoded gradient read-out is to reduce the gradient strength and lower the receiver bandwidth. Decreasing the receiver bandwidth will indeed decrease the noise level, and improve lower signal level detection (proportional to the term kBTB, with kB corresponding to Boltzmann's constant, T corresponding to Temperature in Kelvin, and B is the receiver bandwidth in Hz.) However, this comes at the expense of larger chemical shift artifacts.
Chemical shift artifacts arise as a consequence of the different isochromats associated with different chemical species within the biological sample. In a frequency-encoded k-space read-out, those chemical species which resonate at a slightly higher frequency will appear to be displaced from their actual location in image space towards the direction of increasing frequency. If the spatial frequency encoding gradient is shallow, the apparent displacement can be quite large.
As such, to minimize chemical shift artifacts, the gradient slope is typically made as steep as possible to minimize the apparent shift to within a narrow range (i.e. within 1 or two pixels in the image domain). However, this then requires a larger receiver bandwidth to accommodate the larger frequency range. This in turn increases the overall noise floor at a level proportional to the receive bandwidth.
The conclusion is that frequency readouts generally force a trade-off between gradient strength, noise level, and chemical shift artifacts.
A common technique for noise reduction in signal acquisition is through repeat sampling of a signal and subsequent combination of the multiple measurements. For linear noise sources, such as Gaussian noise, this technique improves SNR through cancellation of the random noise on the signal, the cancellation effect increasing with the number of samples, N.
Noise reduction by this cancellation technique works for static subjects. However, motion-induced blurring is a non-linear effect, so signal combining for which the individual measurements have shifted through large spatial phase angles (relative to the textural/structural wavelengths under study) does not lead to an improved SNR. A fairly standard technique to correct for motion is to look at the MR intensity data in real space and reregister successive traces/images to each other to maximize overlap. It is assumed that, as with the reduction in white noise, linear combination of these reregistered signals will result in reduction of the blurring caused by the motion. However, this only works if the SNR on each individual acquisition is high enough. Reregistering low SNR samples results in a high variance in the estimated position. Threshold theory defines that combining reregistered signals with non-linear blurring, when the original signals are below a certain noise threshold, only increases signal error.
The nonlinearity introduced by subject motion increases at higher k-values, since the motion-induced textural phase shift increases with k—i.e. as the size of the structures of interest decrease, the adverse consequence of motion become more acute. This implies that the multiple samples to be combined need to be derived from the same acquisition sequence, acquired in a sufficiently short time span, to ensure there is negligible motion between samples.
The Cramer-Rao Lower Bound provides insight into the number of samples that are required for a lower bound on the residual variance of an estimate, i.e. the SNR vs. number of samples, in Additive White Gaussian Noise (AWGN). For low source SNRs in AWGN, one needs a large number of samples to average in order to obtain a usable SNR. The primary assumption is that multiple acquisitions can be taken, then averaged to achieve the higher SNR. (CRAMER, H.; “Mathematical Methods of Statistics”; Princeton University Press, 1946. RAO, C. R., “Information and the accuracy attainable in the estimation of statistical parameters”; Bulletin of the Calcutta Mathematical Society 37, 1945.)
Referring to the drawings, the graph in FIG. 1 comparing output SNR shown in trace 102 with number of samples required shown in trace 104 demonstrates that, for high input SNRs, a single sample is sufficient to yield a low noise measure. For lower SNRs, multiple samples are required to “average out” the noise contribution. The ability to combine the samples explicitly assumes that the underlying signal of interest is relatively constant during the multiple sample acquisition process (i.e. the only component which changes is the noise).
The graph in FIG. 2 is a simulation with signal model, trace 202, providing an input SNR, trace 204, showing number of samples of k-value, trace 206, needed to yield a SNR of 20 dB as a function of location in k-space, given an input noise level of 3 mV rms. Since spectral energy density is generally proportional to k−1, to maintain adequate SNR a larger number of input samples is required at higher spatial frequencies (higher k-values). The noise level for the simulation is adjusted for  ̃10 dB SNR at k=2 cycles/mm (λ=500 μm).
As pointed out above, this type of averaging is possible for purely static samples with no displacement or deformation of the targeted tissue occurring over the temporal span of data acquisition. However, for in vivo applications, natural motion occurs even if the patient is compliant. As the texture spacing of interest decreases, the adverse consequences of motion become more acute. More to the point, this type of averaging is based upon the assumption that the underlying signal is the same across acquisitions, and that only the zero-mean, complex-valued, additive, white Gaussian noise (CAWGN) changes. If the signal itself changes, the result will be an average, not only of the noise, but also of the N different versions of the underlying signal, which really doesn't improve SNR.
Using low SNR samples to estimate and correct for motion will result in a high variance of the estimated position. This in turn yields a large variance in the “corrected” acquisitions and does not yield the anticipated increase in SNR when these acquisitions are averaged. This implies that the multiple samples need to be derived from the same acquisition sequence, where motion between samples is extremely small. This is enabled by the method claimed herein.
The issue becomes more acute with shorter structural wavelengths. Consider two acquisitions, noise-free for the moment, one of which has been displaced by an amount d. For a given k-value, an attempt to average them produces:

Y(2πk):=S(2πk)[1+e−j2πkd]/2  (0.1)

Where S(2πk) is the complex-valued signal, and Y(2πk) represents the average of the two acquisitions.

This can be expressed as:

Y(2πk):=S(2πk)e−jπkd cos(πkd)  (0.2)

Which shows both a magnitude attenuation and phase shift, due to the displacement d. Limiting the magnitude attenuation to a floor value a, where 0<a<1, limits d to:





  
    
      
        
          
          ld
          
        
        ≤
        
          
            
              cos
              
                -
                1
              
            
            ⁡
            
              (
              a
              )
            
          
          
            π
            ⁢
            
              
            
            ⁢
            k
          
        
      
    
    
      
        (
        0.3
        )
      
    
  




This shows that, for a given magnitude error, the allowable displacement decreases with increasing values of k. This is because, the smaller the textural spacing of interest, the less motion can be tolerated over the course of data acquisition.
To deal with this problem, an alternate approach is taken in the method claimed herein, which is to dispense with the frequency-encoded readout and to sample specific k-space points, acquiring one or multiple measurements at each k-value of interest at a single spatial location and orientation at a time.
Within a given acquisition in standard MR practice, there are M samples which are acquired of the echo. Instead of acquiring a sample at each k-value, N≦M of those samples could be used for estimation of the (complex-valued) underlying signal value at a specific k-value. Multiple samples within an acquisition can be combined with much less concern of movement than across acquisitions because they are much closer in time.
If the entire echo is used to measure one k-value, the receive bandwidth can be adjusted so as to pass the most abundant resonant peaks in the underlying NMR spectrum, and attenuate frequencies above them.
Taking a straight MRS spectrum (no structural phase encodes), would yield a spectrum consisting primarily of peaks corresponding to H2O (with a chemical shift of δ=4.7 ppm), as well as Carbon-Hydrogen bonds which occur in fat (e.g. CH3, CH2, CH═CH, etc.), each with a different chemical shift ranging from 0.9-5.7 ppm, with the most abundant resonance coming from CH2 in the aliphatic chain which occurs at δ=1.3 ppm.
Assuming use of a 3 T machine, since the Gyromagnetic ratio of Hydrogen is γ=42.576 MHz/T, the chemical shift values are in the range of 166 Hz (for CH2) to 600.3 Hz (for H2O). As long as a (single sided) receiver bandwidth in excess of 600.3 Hz is used, the H2O peak will pass. Assuming baseband sampling, this implies a sampling rate >1.2 kHz (note, if complex base-band sampling is used, this could theoretically be reduced by about 1⁄2.) The point here is that a narrow bandwidth can be used by this method, and sample rates on the order of 800 μs. Noise on the signal is thereby reduced and multiple repeats of the k-value acquisition data are acquired in milliseconds, thereby making the acquired data immune to patient motion. For comparison, a single imaging acquisition is made with a TE of  ̃30 ms, and TR on the order of 500 ms-2000 ms. To acquire the repeats necessary for signal averaging can take minutes—a temporal range wherein respiratory, cardiac, and twitching motion limits resolution through motion-induced blurring. The claimed method enables acquisition of values in regions of k-space which have very low signal levels, such as would be found for higher k-values (shorter textural-wavelengths)—the fine texture range that has hitherto remained elusive.
To maximize the signal, the non-zero frequencies of abundance are selected. In general, this does not correspond to a mere averaging of all of the samples acquired. Instead it is akin to a matched filter which is “tuned” to the frequency of interest, corresponding to the specific chemical species of interest.
As a side note, the full NMR spectrum may be extracted (without any phase encoding gradients: just volume selection) to obtain a baseline of the underlying signal strength (and associated frequencies), which in turn will be spatially modulated, providing insight into textural wavelengths through knowledge of the chemical species expected in the textural elements under study.
The isochromats of interest can be extracted by acquiring N samples of the echo, then taking the Fourier transform. Since the echo is being played out with no gradient, the strength of the resulting signal at the Isochromat of interest will correspond to the (complex-valued) k-value coefficient of interest.
Given the goal is to extract the relative magnitude of textural wavelengths, just the magnitude vs. textural wavelength measurement is the required information. However, in order to extract sufficient signal strength and differentiate it from the underlying noise floor, the complex phasor values must be preserved until the end.
The relationship between the noise floor, the signal strength (at a specific isochromat where there is an abundance of chemical species), the number of samples required, and the max tolerated error can be approximated as




  
    
      
        N
        ≥
        
          
            σ
            2
          
          
            
              
                
                A
                
              
              2
            
            ⁢
            
              ɛ
              2
            
          
        
      
    
    
      
        (
        0.3
        )
      
    
  




Where σ2 represents the noise variance, |A|2 represents the squared magnitude of the isochromat(s) of interest, and 0<ε<1 represents the allowable error of the estimate. Further assuming that the noise is mostly sourced from the biological sample, this can be further approximated as:





  
    
      
        N
        ≥
        
          
            
              
                N
                eff
              
              ·
              
                k
                B
              
            
            ⁢
            TB
          
          
            
              
                
                A
                
              
              2
            
            ⁢
            
              ɛ
              2
            
          
        
      
    
    
      
        (
        0.4
        )
      
    
  




Where NFeff is the effective noise figure of the receiver, kB is Boltzmann's constant, T is the temperature in Kelvin of the biological sample, and B is the receiver bandwidth. In this case, N can be used as a guide to the number of samples that need to be acquired within a given acquisition in order to create a reasonable estimate.

If the number of samples required exceeds the number available in one acquisition, combination of measurements from a single acquisition may be needed to maximize the signal, prior to spatial reregistration between acquisitions. A reasonable estimate and displacement correction between the two or more acquisition sets is needed. Combination of measurements at a single k-value from a single TR block can now be used to boost the SNR such that reregistration between successive TRs has a much greater chance of success.
While the entire set of samples acquired in an echo or entire TR could be allocated to the estimate of one coefficient in k-space, if acceptable values can be estimated using fewer than the maximum number of echo samples, it opens up the possibility of being able to acquire more than one coefficient in k-space within a specific echo or TR.
FIG. 3 shows an example timing diagram for a pulse sequence for data acquisition using the method claimed herein. RF pulses included in trace 302 are employed to excite selected volumes of the tissue under investigation, as in typical MR imaging. A first RF pulse, 304, is transmitted coincidentally with a gradient pulse 308 on the a first magnetic field gradient, represented in trace 306. This excites a single slice, or slab, of tissue the positioning of which is dependent on the orientation and magnitude of the first gradient, and the frequencies contained in the RF pulse. The negative gradient pulse, pulse 310, rephases the excitation within the defined thickness of the slice or slab.
A second RF pulse 312, at twice the magnitude of first RF pulse 304, is transmitted coincidentally with gradient pulse 316, on a second gradient, represented in trace 314, exciting a slice-selective refocus of spins, this second tissue slice intersecting with the first. (As this second RF pulse 312 tips the net magnetic vector to antiparallel to B0, it results in inversion of spins and subsequent refocusing, thus leading to a signal echo at a time after the 180 degree RF pulse equivalent to the time between the 90° and 180° RF pulses.) An initial higher value gradient pulse, 318, at the start of gradient pulse 316 is a crusher, or “spoiler” gradient, designed to induce a large phase wrap across the tissue volume. A similar gradient pulse, 322, at the trailing end of pulse 316, as it comes after the 180 degree RF inversion pulse, unwinds this phase wrap. In this way, any excitation that is not present prior to the 180 degree RF pulse, such as excitations from imperfections in the 180 pulse itself, will not have this pre-encode so will not be refocused by the second crusher, hence will not contribute to the signal. In summary, the second RF pulse, in combination with the applied second gradient, provides slice selective refocusing of the signal in a region defined by the intersection of the first slice and the second slice set by this second gradient.
An encoding gradient pulse 326, on trace 314, sets an initial phase wrap, hence k-value encode, along the direction of gradient pulse 326. In general, the k-value encode can be oriented in any direction, by vector combination of the machine gradients but for ease of visualization is represented as on the second gradient.
A refocusing third RF pulse 328, applied in combination with gradient pulse 332 on a third gradient, represented by trace 330, defines a third intersecting slice selective refocus to define the VOL Gradient pulse 332 again employs crusher gradients.
The negative prephasing gradient pulse 326 winds up phase such that, at the signal echo following the second 180° RF pulse, signal acquisition starts at high k-value, which may then be subsequently decremented (or incremented or varied in orientation) for further acquisitions, as will be described below. As energy density in the signal is generally proportional to k−1, this method ensures k-values with lower SNR are acquired first, before T2 effects have caused much overall signal reduction.
With all gradients off, a receive gate 333 is opened to receive the RF signal, which is shown in FIG. 3 as pulse 334 on trace 336. The RF signal in trace 336 is a representation showing only the signal present in the receive gate window without showing the actual details of the RF signal outside the window. Sampling occurs as represented by trace 338 beginning with the initial k-value, 340a, seen on trace 324. Note that, at the scale of the drawing, the sampling rate is high enough that the individual triggers of the analog to digital converter (A/D) have merged together in trace 338. (The expanded time scale in FIG. 4 described below shows the individual A/D triggers.)
In regions of k-space where the corresponding coefficients are sufficiently large that they can be well-estimated using a small subset of the samples of one echo, acquisition of another k-value, obtained by applying a gradient pulse 342a shown on trace 314, to select a new k-value, during the time the echo is being recorded, is accomplished. After a suitable settling time, another set of samples of the echo (now derived from the new k-value coefficient) can be collected. This process can be repeated, acquiring multiple samples at each of a select set of k-values within one TR. A plurality of samples are taken at the initial k-value 340a. A k-value selection gradient pulse 342a is then applied and the resultant k-value 340b is sampled. (Though shown in the figure as a negative pulse on the second gradient, decrementing the k-value, in practice this pulse and subsequent k-value gradient pulses can be designed through any vector combination of gradients to select any k-value or orientation.) Similarly, the k-value selection gradient pulse 342b, selects a third k-value 340c which is sampled by the A/D. Each gradient pulse changes the phase wrap, selecting a new k-value. Application of a k-value selection gradient pulse (342c-342f) followed by multiple sampling of the resultant k-value coefficient is repeated as many times as desired. While data is being acquired throughout, the samples of interest are acquired when all gradients are off. The gradient orientations for slice and k-value select may be coincident with the machine gradients, which are aligned to lie coincident or orthogonal to the B0 field. Alternatively, the acquisition directions and k-value encodes may be selected using gradients that are a vector combination of all three machine-gradient axes.
In the circumstance where it is desired to measure a low SNR k-value the prewinding encoding gradient pulse can be set such that the first k-value to be measured is the desired low SNR k-value. Alternatively, the prewinding gradient pulse can be set to zero so that the first k-value measured is k0. A measurement of k0 may be desired for the purpose of determining the systems receiver sensitivity to the particular VOI, determining the relative prevalence of isochromats (e.g., water vs. lipids) irrespective of texture in the VOI, or for the purpose of establishing a reference value for normalization of the other k-values measured in a VOI or for comparison with k-values from other VOL Furthermore a strategy for gathering a specified set of k-values for a VOI may include measuring the low SNR k-values (typically the higher k-values) in a first set of multiple TR and then measuring k0 and other higher SNR k-values in other TRs while remaining in the same VOI.
As is shown diagrammatically in FIG. 3, the signal magnitude reaches a maximum at the time of the spin echo. It is also shown diagrammatically that the signal magnitude is varying throughout the acquisition of the multiple RF measurements of a k-value and more so between successive blocks of measurements of k-values. Alignment in time of the measurement of the low SNR k-values with the highest echo signal enhances the SNR of the k-value measurement, alternatively alignment of higher SNR k-values with lower echo signal magnitude allows gathering additional useful k-value acquisitions during the echo.
FIG. 4 shows a close-up of the pulse sequence of FIG. 3 during the initial portion of the RF sampling window 338 between 7.25 and 8.00 msec. Multiple samples of the same k-value, taken in rapid succession with all gradients off, provide the input for signal averaging to reduce AWGN when SNR is low. In a first block 344a of the sampling window 338, multiple samples 346a are taken of the first k-value 340a. During application of the k-value selection gradient pulse 342a, transition samples 348a are taken. When the k-value selection gradient is switched off, multiple samples 346b are taken at the second k-value 340b. Application of k-value selection gradient pulse 342b then occurs with associated transition samples 348b, and subsequent acquisition of samples 346c of the third k-value 340c after the gradients are switched off. The underlying signal is minimally impacted by motion due to the very short time window used to acquire data at each given k-value. Since the data is acquired with gradients off, there is no issue with chemical shift and the effective T2* is longer, boosting the signal value.
The sampled values of the echo, acquired while the k-value selection gradient pulse is ramped up, held steady, and then ramped down to zero, will necessarily be influenced by the applied gradient. These transition samples may provide other interesting information, but are not used in the consideration of a straight measurement of the k-value coefficient; only those samples which are recorded when there is no gradient currently active are used for this.
A consistent number of samples at each k-value can be acquired, or an alternative sequence may be employed where, as k-values decrease, hence increasing in signal amplitude, fewer samples are acquired. A pulse sequence designed for this type of acquisition is illustrated in FIG. 5. Multiple samples of each k-value targeted in the acquisition are acquired in rapid succession, with all gradients off. These repeats provide the input for signal averaging in low SNR signals. As with the pulse sequence, depicted in FIGS. 3 and 4, the underlying signal is minimally impacted by motion due to the very short time window in which data is acquired for a given k-value.
Samples within the portions of the sample window 344a-344g outlined on FIG. 5 correspond to the number of samples acquired for a given k-value 340a-340g each induced by an unwinding pulse 342a-342f of the k-value selection gradient. Nk, the number of samples associated with a given k-value, can be selected based upon expected SNR, tissue contrast, contrast to noise, pathology, texture size, and/or texture bandwidth. For the example in FIG. 5 it can be seen that a decreasing number of samples are taken for progressively smaller k-values (larger textural features). This is because, as previously discussed, to first order signal amplitude increases with decreasing k-value—energy density is generally proportional to k−1. For this same reason, larger k-values are acquired first in this scheme, when T2 effects are least, the longer wavelength, higher signal strength, k-values being recorded later in the acquisition.
Refocusing the echo, and/or a new TR can be used to build up a library of k-space samples. Acquisition of multiple k-values within one TR can be facilitated by application of multiple refocusing gradients and/or RF pulses, to increase the time over which the additional k-values can be sampled within a TR. These later echoes would presumably be used to acquire the coefficients of the lower k-values in the selected set, as their energy density in the continuum of values is generally higher so the effect of T2 decay on overall signal magnitude will not affect them as severely as it would the higher k-values. In this way a larger portion of the required k-space filling can be accomplished over fewer TRs, allowing more rapid data acquisition, minimizing the need for repositioning the VOI.
FIG. 8 shows an extension of the basic sequence of the method claimed herein, using spin-echo refocusing to extend the record time for the TR. Application of a refocusing RF pulse 802 with an associated gradient pulse 804 results in slice-selective refocusing. After an appropriate settling time, a second sampling window 806 is opened by the receive gate 808. Multiple k-value selection gradient pulses 810 are applied to increment the selected k-value and, after switching off each successive gradient pulse, multiple samples of the selected k-value are acquired in the sampling window. A second slice-selective refocusing RF pulse 812 with associated gradient pulse 814 again inverts the spins and, after application of each in the multiplicity of k-value selection gradient pulses 820, data is acquired in a third sampling window 816, opened by the receive gate 818. As shown in the drawing, an increasing number of k-values may be sampled with each refocusing. Refocusing can be repeated until the decrease in signal level from T2 and other effects makes further signal acquisition ineffective. Another method to extend the record time by exciting multiple signal echoes, is to use one, or a series of, gradient recalled echoes (GRE). GRE are different from the SE in that they cannot refocus the effects of stationary inhomogeneities, so T2* effects limit the number of repeats.
In addition to the tissue contrast available, the k-values associated with particular pathology will be part of the determination of the number of samples needed for signal averaging, Nk. In liver fibrosis, the wavelength of pertinent textures is in the range of 400 microns, i.e. a k-value of 2.5 cycles/mm. This is similar to the textural spacings seen in fibrotic development in many other diseases, such as cardiac fibrosis. The spacing of elements in trabecular bone varies a lot, but the minimum spacing of interest is the width of trabecular elements, which are approximately 80 microns, setting a maximum k-value of 12.5 cycles/mm. In neuropathology, many of the textures of interest are very fine, on the scale of 50 microns, equivalent to a k-value of 20 cycles/mm
Each pathology will dictate what exactly is needed as quantitative data, i.e. what part of the continuum of k-values needs to be monitored, and with what resolution and sensitivity. In some pathologies, short (long) wavelength features increase at the expense of long (short) wavelength features (e.g. liver fibrosis). In other pathologies, an amplitude decrease and broadening of short wavelength features indicates disease progression—e.g. degradation of the ordered formation of cortical neuronal minicolumns (approximately 80-micron spacing) with advancing dementia. In bone, with increasing age, first the highest k-value features disappear in the structural spectrum. Next the major structural peaks shift slowly towards lower k-values with advancing osteoporosis, the pace of this shift accelerating as an increasing percentage of trabecular elements thin to the point that they break.
The signal level obtainable will depend on anatomy to some extent. For instance, though the resolution needed is highest in brain, the proximity of the cortex to the surface of the head ensures that use of a surface coil will provide significant signal boost for cortical structures. Lower resolution is required in liver, as the structures of interest are on the order of several hundred, rather than tens of microns. But, the organ is deeper (further from the coil) reducing the measured signal. Using the in-table coil for spine data acquisition yields modest signal level and good stabilization. Also, bone is a high contrast target, so the SNR requirement is not as stringent. For all these reasons, the exact number of repeats needed for averaging depends on more than the k-value range targeted.
FIG. 6 shows a simulation demonstrating that the ability provided by the claimed method to acquire many repeats of a targeted k-value within a single TR enables robust signal averaging to boost SNR. Assuming a subject displacement rate (which has in practice been measured clinically over the course of several scans) of 30 μm/sec, and a sampling rate=33.3 kHz (ΔTsample=30 μs), 90 repeat samples for averaging can be taken rapidly enough that, even up to a k-value of 20 cycles/mm (texture wavelength=50 μm), the acquisition remains immune to motion effects.
FIG. 7 shows that, for comparison, using the conventional approach of acquisition of a spatially encoded echo, even assuming a relatively fast gradient refocus sequence, which would provide a sampling rate of about 67 Hz (ΔTsample=15 msec), subject motion over the time needed for 90 repeats would severely degrade the signal, and any ability to improve SNR by signal averaging. The situation is actually worse due to the fact that to acquire 90 repeats using conventional spatially-encoded echoes would require several TRs, making the acquisition time significantly longer, and the signal degradation due to motion much more severe. With the exception of the very lowest k-values, the potential SNR gain due to multiple sample combination has been nullified by the effects of motion.
By acquisition of a large enough selected range of k-values, construction of a structural profile in one or more dimensions becomes a possibility. As discussed above, refocusing echoes within a single TR, or multiple TRs, can be used to build up a library of k-space samples. Phase coherence might not be maintained between different k-values if they are acquired in TRs separated temporally such that displacement has occurred between them. If our primary interest is in the relative strength of signals at particular k-values, this is not a problem. If creation of a profile or an image from this library of values is desired, the necessary post processing will have as input the high SNR measure obtained within each TR using the method claimed herein. These measures can then provide robust input for any required reregistration between echoes or TRs towards constructing a profile. As an example selection and measurement in a first TR of a set of selected k-values may be accomplished with at least one having a low k-value. In a subsequent TR, selection of the same set of k-values will allow re-registration of the data between the two TRs since even if significant motion has occurred the phase change in the low k-value phase shift will be less than for higher k-value textures and may be correlated between the two TRs. Basically, the higher the k-value, the greater the phase shift due to subject motion. Acquiring signal from successive encodes with a large difference in k-value enables a better estimate of phase shift by careful comparison of the apparent phase shift for each.
This is very similar to x-ray diffraction, wherein the magnitude-only information (no phase) obtained presents the challenge of determining a best estimate of the corresponding structural profile based on this magnitude-only information. Algorithms exist towards solving the problem, the chance of success depending on the range of k-value coefficients obtained, the SNR of each averaged coefficient, and the width of values contained in a nominally single-valued acquisition of k-value. The chance of success in this effort is greatly increased using the claimed method due to its immunity to subject motion.
The ability to reconstruct a profile from k-value data depends somewhat on the spectral broadness of each single-k-value acquisition. While this is influenced by the VOI (Volume of Interest) size and shape, it is also influenced by k-value and pathology, as degradation of tissue often tends to lead to more textural randomness within tissue.
Selection of the VOI—shape, dimensions, orientation, and positioning within an organ/anatomy affects the data measured and its interpretation. The VOI shape can be chosen to maximize usefulness of the acquired data. Data can be acquired in different directions, and at different textural wavelengths (k-values) within a VOI enabling assessment of textural anisotropy. Texture can be sampled in multiple VOIs, either interleaved within a single TR, or in successive TRs, towards assessment of pathology variation across an organ. Standard interleaving processes for the VOI may be used within a TR to provide additional data by applying additional encoding pulses on vector combination gradients and associated k-value selection gradient pulses for k-values in the interleaved VOL As previously described, additional excitation RF pulses with associated slice selection gradients may be repeated within the same TR by exciting a volume of interest with a gradient set in each repeat having at least a first gradient with an alternative orientation from the first gradient pulse 308 applied initially in the TR, to define an additional VOI for excitation in new tissue, that does not overlap any previous VOI in the TR (fourth, fifth and sixth gradients in a first repeat and succeeding incremental gradients in subsequent repeats). This response can be mapped, or the several measures taken and averaged, whatever is appropriate for the targeted pathology. This is similar to the multi-positioning of tissue biopsy. However, in the case of tissue biopsy, the number of repeats is limited due to the highly invasive nature of the technique. The minimum number of structural oscillations to be sampled at a specific k-value dictates a minimum VOI dimension in the direction of sampling—the length required varying inversely with targeted k-value.
To ensure adequate sampling of structure when targeting a range of k-values, the VOI dimension in the sampled direction can be held constant for all k-values in the targeted range, with the result that the number of structural oscillations sampled will vary with k-value. This is a simple solution, requiring the sampling dimension be set by the lowest k-value (longest wavelength structure). Using this approach, the sampling dimension of the VOI is larger than required for the highest k-value in the range, thus providing less localization within the tissue than would be otherwise possible.
Alternatively, data at widely differing k-values can be acquired in successive TRs, using changing VOIs tailored to the specific k-value targeted. Or, the dimensions of the VOI can be selected such that acquisition in different directions within the VOI will be tailored to sampling in a specific textural frequency (k-value) range.
Similarly, the VOI may be held constant and the vector combination gradient for the encoding and k-selection pulses may be altered from TR to TR for assessing feature size.
In some instances, it is desired to localize tightly in the spatial domain, to broaden the localization in k-space. By defining a non-cubic acquisition volume, it would be possible to acquire data from differing k-values along the different (orthogonal or other) directions within the VOI, within one TR. The elliptical cross-section VOI 902 in FIG. 10 is one such possibility. Acquisition along any radial direction, as well as along the axis of the shape, would be possible within one TR.
Additionally, one could use the flexibility of the method claimed herein to sample k-space in a linear or in a curved trajectory. For example, texture could be sampled along radial lines, or along an arc or a spiral, to extract information of textural sizes along different spatial directions. These methods can be used to determine the anisotropy of texture, or the sensitivity to alignment in structures that are semi-crystalline, such as cortical neuronal columns, or to more rapidly build up a library of k-values within a targeted extent of tissue in an organ.
During one TR (i.e., one 90 degree excitation) k-value encodes can be applied in multiple directions by changing the applied vector combination gradients for encoding and k-selection pulsing. The exact form of the VOI and sampling direction can be used to yield much textural information. For instance, the organization of cortical neuron fiber bundles is semi-crystalline, as the bundles in healthy tissue form in columns. Because of this, the measure of textural spacing perpendicular to the bundles is very sensitive to orientation. When the orientation is exactly normal to the columns, a very sharp signal maximum is expected, the signal falling off rapidly as the orientation varies on in either rotational direction away from this maximum. One way to measure the spacing and organizational integrity (a marker of pathology) would be to “rock” the acquisition axis around this maximum looking for a resonance in signal intensity. This approach of looking for “textural resonances” by looking for signal maxima can be applied in any tissue region. As pathology degrades the organizational integrity, the sharpness of this peak will degrade and the signal maximum will be reduced.
Similarly, the randomness of the spacing in certain textures can be assessed by varying the length of tissue sampled in a specific, or in multiple directions, with subsequent change in acquisition length. The selected value for that length can be varied over multiple TRs to test the sensitivity of the measured coefficient to this parameter.
The VOI can be selectively excited by a number of methods, for instance intersecting slice-selective refocusing, selective excitation using phased-array transmit in combination with appropriate gradients, adiabatic pulse excitation to scramble signal from the tissue outside the region of interest, as examples. Parameter selection for the various methods can be done with SNR optimization in mind. For instance, the VOI generated by a slice selective excitation and two additional mutually-orthogonal slice selective refocusing pulses, as by VOI 904 in FIG. 9. Through careful RF pulse design, the shape of the VOI can be designed so that the edges are smooth and more approximate a windowing function, as shown in FIG. 10. These windowing functions provide the volume selection without adverse impact on the spatial frequencies. Recall, in Fourier theory, each spectral line is smeared by the convolution of a Fourier transform of the window function. It is desirable to minimize this smearing of the underlying spectrum, as it decreases the energy spectral density, and adversely impacts SNR.
Importantly, as has been discussed previously, the VOI can be moved from place to place within an organ or anatomy under study to measure the variation of texture/pathology. This response can be mapped, or the several measures taken and averaged, as appropriate for the targeted pathology. This is similar to the multi-positioning of tissue biopsy. However, in the case of tissue biopsy, the number of repeats is limited due to the highly invasive nature of the technique.
Different diseases and conditions affect tissue in different ways. Generally, pathology advancement entails: 1) a loss of energy density in specific regions of k-space, and/or 2) a shift in textural energy density from one part of k-space to another, both effects being accompanied by 3) changes in the width of existing peaks in the continuum of textural k-values. Using trabecular structure as an example—with decreasing bone health the average separation of trabecular elements widens (texture shifts to lower k-values) and becomes more amorphous (broader peaks in k-space), while in parallel the structural elements thin (a shift to higher k-values in a different part of the spectrum). Other tissues/organs are affected by diseases that have their own individual signatures in tissue texture.
Using the method claimed herein, k-space is probed to reveal texture in such a way as to eliminate the loss of signal resolution that results from subject motion blurring. Instead of measuring the large continuum of k-values needed to create an image, the focus here is on acquisition of a select few k-values per TR, with sufficient repeats of each to yield high SNR. Each of the individual acquisitions is centered on a single k-value. While the spatial encode is, to first order, a single spatial frequency sinusoidal encode, there are a number of factors which have the effect of broadening the spatial frequency selectivity of the k-value measurement. One significant factor affecting the broadness, or bandwidth, of the k-value measurement is the length of the sampled tissue region. A longer sampling length encompasses more textural wavelengths along the sampled direction, which has the effect of narrowing the bandwidth of the k-value measurement. (This is the inverse relationship between extent of a measurement in real and in k-space.) Hence, an aspect of the claimed method is the ability to set the bandwidth of the k-value measurements by appropriate selection of the sampling length determined by the VOI dimensions or determined by the acquisition dimensions. Using this method, the bandwidth of the measurement can be set according to the desired k-space resolution appropriate to the tissue being evaluated. (Need both high k-values for good texture resolution, and high resolution in k-space for sensitive monitoring of pathology-induced changes.) For highly ordered structures one could choose a set of narrow bandwidth measurements distributed over the expected range of texture wavelengths, whereas in a more randomly ordered structure, such as the development of fibrotic texture in liver disease, one could choose to use a single, or a few, broadband k-value measurements to monitor development of the fibrotic texture.
A measure of both the relative intensities of the various textural k-values present in a tissue and the broadness of the peaks along the continuum of textural k-values present within the texture under study is needed. As such, data acquisition can be designed to probe specific region(s) of k-space, with parameter selection that will enable measurement of the relative width of peaks arising from the underlying tissue, rather than that resulting from experiment parameters. It is necessary to recognize the interaction of the two components, and design experiments to yield the best measure of pathology-induced tissue changes.
It is desirable to obtain a good measure of texture by acquiring multiple measures of signal amplitude at specific k-values close in time before motion blurs the data, taking repeat measurements in minimal time to allow best inter-measure correlation for averaging. an alternative to acquiring many repeat measures at one point in 3D k-space, is to acquire data with a gradient on, such that the k-value is changing continuously across the acquisition, the extent in k-space being determined by the height of the gradient and its pulse width. In addition to varying the magnitude of the k-vector, its direction over the course of data acquisition can also be varied. Combination of direction and magnitude changes across an acquisition result in a curvilinear trajectory through k-space. If this deviation is small enough that the k-values remain correlated to some extent, they can be combined more effectively to increase SNR than if they were simply averaged. Gradient on acquisition can therefore be used to intentionally vary the direction and magnitude of the k-vector, for the purpose of smoothing signal speckle—which manifests as a time varying signal over the data acquisition, resulting from interference of the individual spin signals' varying phases and amplitudes. The selected variation in k-value direction and magnitude across the acquisition is chosen to provide sufficient combined measures to get an estimation of the representative power within a neighborhood of k-space.
Varying the k-value to reduce speckle can be accomplished within a single, or multiple, echoes. For a sphere in k-space, defined by the magnitude of the k-value under study, the k-value can be varied by keeping the magnitude of k-constant but sweeping the vector over the surface of the sphere, or the same angular orientation may be maintained, and the magnitude of k varied, or both can be varied simultaneously
For the purpose of reducing speckle effects, these variations would usually be small enough deviation from either the k-magnitude or direction that there is a meaningful correlation between the measurements for the particular tissue under investigation.
The major components of the spatial frequency will be the same in all those measurements (they are correlated) unless the measured tissue is a highly crystalline texture. But the normal diffraction pattern for a micro-crystalline or amorphous structure has a lot of speckle. Consequently, by sampling a number of points in the same region of k-space they can be combined in various ways, selected to provide optimal smoothing to reduce the speckle-pattern. A better and more robust measure, from averaging out the fluctuations, is the result.
A number of approaches to “dither” k-value to reduce speckle or to tailor width in k-space may be employed. A first approach employs constant k-magnitude plus sweeping through a range of angles by keeping gradients on during acquisition and combining the measures using correlative information to eliminate speckle. Alternatively, the same direction in k-space may be maintained but the magnitude varied by leaving gradient on during acquisition and combining the measures using expected correlation. As yet another alternative, both magnitude and direction may be varied simultaneously or over an acquisition series, essentially performing the other two alternatives simultaneously to both reduce noise and provide a better assessment of the representative k magnitude in a structure in a “small” region around a specific k-value, i.e., to reduce speckle.)
For combining the measures at different k-magnitudes, for noise reduction averaging, there is a phase shift from one radius (magnitude) to the next from the gradient wind-up. Rephasing may be accomplished before averaging.
Combining the different magnitude measures in an amorphous structure is more well-known than combining different angled measures. Now in addition to the scheme of reducing thermal noise by rapid sampling the fluctuations due to speckle (which though real signal confounds good assessment of the spatial frequency) may be reduced.
A dynamic k-space acquisition is therefore employed. The acquisition mode is dynamically chosen based upon the Signal to Noise Ratio (SNR) of the signal at various k-space locations. The gradient, applied during signal acquisition, post-acquisition receive bandwidth, and estimation algorithms used are dynamically adjusted based upon the expected SNR values in k-space to optimize acquisition time and post-processed SNR. In regions of high SNR, a single sample at a given k value may be a sufficient estimator. This requires a relatively wide receive bandwidth to accommodate the relatively rapid signal variations in the receive chain as k is changed rapidly (due to the large gradient).
In regions of moderate-to-low SNR, the gradient magnitude is decreased so that, subsequent samples, while not taken at identical k values, are correlated, which in turn can be used to improve estimates of the underlying signal values within that range of k-values.
Correlation may be introduced in k-space due to selected windowing in profile space. To enable combination of sequential samples from the ADC so as to improve SNR, correlation among successive samples will be increased by proper choice of windowing in profile space, a shorter window driving a greater correlation distance across sequential values in k-space and a longer window resulting in lower sample to sample correlation. Inducing correlation of neighboring points in k-space by windowing in profile-space is a mathematical tool that can, in many cases, help to measure the underlying texture in a low SNR environment. Basically, windowing blurs the data so that the k-value power spectrum is smeared out through k-space, so that sequential measures can be averaged/combined more easily to increase SNR.
In a very high SNR environment as large a window as possible is used because is a measure of the actual textural power distribution across a range of k-space is desired. The longer the sampled region in real (profile) space the more accurate the measurement when measuring amorphous textures. Reducing the sampled region by windowing to induce correlation in k-space actually obscures the specific desired measurement point.
However, while facilitating measurement, inducing correlation through windowing does blur to a greater or lesser degree the underlying relative power density profile in k-space arising from the underlying texture, which is the target of the measurement. As the sample-to-sample spacing (determined by the analog to digital converter speed and gradient height) in k-space decreases, there will be increased correlation, which can be used in post-acquisition processing to form better estimates. Additionally, the receive bandwidth in these regions can be decreased, which further decreases the noise floor.
In regions of very low SNR, multiple acquisitions of the measured signal level at a specific k value can be taken, with zero (and/or non-zero) gradient during acquisition. The multiple acquisitions can then be optimally combined to provide an estimate for specific k values.
An example Hybrid pulse sequence is shown in FIG. 11. This particular sequence is an example of a rapid acquisition with refocused echoes (RARE) type sequence wherein three different levels of gradient are used for acquisition as illustrated across three separate echoes. The shown pulse sequence for selecting the desired VOI and initial phase wrap for k is as described in FIG. 3 and is numbered consistently in FIG. 11. While this exemplary pulse sequence for establishing the VOI is employed in the various examples disclosed herein, the determination of VOI made be made by any of numerous approaches including, as an example, time varying RF pulses with commensurately time varying gradients applied. Similarly, during or after the pulse sequence employed for determining the VOI, the encoding gradient pulse may be applied for selection of the initial k-value. The data recording starts with acquiring values in regions where |k|>>0, given that the signals are smallest there and should be acquired first. The second echo samples values associated with |k|>0, but whose signal levels are still relatively small and require combination of multiple measures to provide robust SNR. The final echo samples values associated with |k| in the neighborhood near |k| ̃0 where the corresponding signals are largest. Note that this is just one example of how this hybrid approach, using both zero and non-zero gradient in one acquisition, could be used. Different amounts of k-value windup (as determined by the gradient height and pulse duration) can be acquired in one echo rather than in multiple echoes as will be described subsequently. Multiple combinations of the differing k-value windup also can be acquired within one echo. Additionally, while refocusing is disclosed in the drawings as employing an RF pulse, gradient refocusing may also be employed.
Details at expanded scale of the pulse sequence/signal acquisition are shown in FIG. 12. In this sequence, the first echo from RF pulse 328 with gradient pulse 332 uses the a pulse sequence similar to that previously described with regard to FIG. 3 and acquisition of multiple samples with gradient=0 and incrementing of k values with k value selection gradient pulses 1100a, 1100b and 1100c. This samples multiple values of a given location in k-space, which values are then optimally combined. This is appropriate for regions of k-space whose values are very small and therefore have very low SNR. This typically occurs in regions where |k| is large.
The second echo from RF pulse 1102 with gradient pulse 1104 is acquired with a small non-zero gradient 1106 acting as a time dependent phase encode. A small gradient may be defined as a gradient inducing samples in k-space which will be sufficiently closely spaced so that the samples are highly correlated. These samples can then be post-processed by an estimator which takes advantage of the high inter-sample correlation to improve the resulting SNR. Quantitatively an exemplary “small” gradient might be up to 20% of the magnitude of the encoding gradient pulse. As seen in the figure samples of multiple values in a relatively small neighborhood, Δk, in k-space are obtained. The spacing of Δk can be chosen such that, due to the windowing of the VOI, there is high correlation between neighboring samples. The correlation is exploited in the estimation algorithm to generate an optimal estimate of the signal levels across a neighborhood in k-space. This is appropriate for regions of k-space which have low SNR, but whose values, because of the correlation induced by the windowing function, vary slowly across k-space in a small neighborhood.
A third echo from RF pulse 1108 with gradient pulse 1110 is acquired with a relatively larger time dependent phase encode gradient 1112. The higher gradient employed herein creates sequential measurements in k-space which have a lower degree of correlation across the neighborhood of k-space under study. In this case, there is lower inter-sample correlation available for SNR improvement. A higher gradient may be employed for sampling k-space locations whose values have a high SNR to begin with (such as would be seen in lower textural frequency (low k-value) regions). As seen in the figure samples relatively widely spaced across k-space, well outside the inter-sample correlation imposed by the window function, are generated over the entire pulse. This is appropriate for rapid acquisition of values in k-space whose signal levels are relatively high and enable high SNR recording. In this case, a single sample at a given point in k-space provides a high enough signal. However, in such higher SNR regions, the higher gradient may be employed and selected bursts of data samples may be rapidly recorded. Each of these bursts will have substantial inter-sample correlation within the burst and may allow computation of results similar to that described for the lower gradient acquisition discussed above. This may be viewed as transitioning along the K-space line with the higher gradient while sampling blocks of data at closely spaced k values to maintain correlation.
Non-zero gradient acquisition allows sweeping across a curvilinear path in k-space. By a judicious choice of time increment, At, gradient magnitude, G, w(X), and the total number of samples acquired, N, the neighborhood of “high” correlation can be adjusted to be M<N. This in turn would allow estimation of a multiplicity of distinct values within k-space by using a subset of M samples for each output estimation. Textural data from within a tissue region defined by the VOI can be acquired with the non-zero gradient to enable determination of the local distribution of power density of k-values within a neighborhood in k-space. The extent of k-space sampled with a gradient pulse played out during acquisition is determined by the gradient height and the gradient pulse width (pulse duration). The spacing between signal samples in k-space is determined by the gradient height and the sampling rate (limited by the maximum speed of the analog to digital converter). The correlation between sequential samples in k-space is determined then by the spacing between samples, by subject motion, by the window used to bracket the acquisition in physical space, and by the underlying texture.
A useful method for selecting the acquisition parameters is with reference to the degree of correlation needed within a set of values to be combined. The wavelength of a repeating structure (texture) is defined as the inverse of the k-value associated with that texture, λtexture=1/ktexture. To be able to combine a set of values [measurements] to yield improvement in SNR, the underlying textural signals must not be shifted in phase by a significant percentage of λtexture relative to each other. In exemplary embodiments, the phase shift across the set of samples to be combined should be no greater than 80% of 2π.
Resolution in MR imaging is limited by subject motion during image acquisition. This limitation can be very severe with non-compliant patients. In addition to patient motion/compliance, the resolution achievable in MR imaging which is exemplary art comparable to the present invention, depends on several factors, such as tissue contrast, organ, coil type, proximity to coil. Robust imaging of structures below about 5 mm in extent is problematic, and anything below about 1 mm is outside the realm of routine clinical imaging. This is a clear shortcoming as many tissue textures in the range of about 5 mm down to 10 μm develop and change in response to pathology development, hence measurement of these textures can provide much diagnostic information—these tissue changes are most often the first harbinger of disease. It is this textural wavelength range, from about 5 mm down to 10 μm that is targeted with the presently disclosed method.
To measure tissue texture, the range of wavelengths in real space which can be resolved, i.e. the wavelengths of the textures pertinent to the particular pathology, are in the range of several mm down to microns. This is the range made inaccessible (blurred) in imaging due to patient motion. As k is defined as 1/wavelength, a range of k-values from about 0.2 mm−1 to 100 mm−1 is employed in exemplary embodiments to define the textures of interest. This brackets the region of k-space of interest, and defines the gradient height and duration of the encoding gradient pulse to induce phase wrap to create a spatial encode for the specific k-value and orientation as well as for the non-zero gradients applied for measurement of the neighborhood around the initially selected k-value. The method of the embodiments herein for sample acquisition and post processing may all be conducted in k-space. The only localization in real space is the positioning of the VOL Just enough of the neighborhood around a point in real-space is sampled to measure texture—i.e. to determine the power distribution within a neighborhood in k-space around the selected point in real space.
The exact range needed varies with the targeted pathology. For example:
Osteoporotic development in bone microarchitecture. As examples, the variation in average trabecular spacing (TbSp) from healthy to osteoporotic bone brackets a wavelength range of about 0.3 mm to 3 mm; the equivalent range of k-values is 0.34 mm−1 to 3.4 mm−1 With fibrotic liver disease monitoring change in liver tissue texture from the healthy collagen-highlighted vessel-to-vessel spacings to the diseased state in which the lobule-to-lobule spacing becomes the prominent tissue texture. Vessel-to-vessel range is 0.4 mm to 1.5 mm translating to k-values of 0.67 mm−1 to 2.5 mm−1 while lobule-to-lobule spacing of approximately 1 mm to 4 mm, translating to k-values from 0.25 mm−1 to 1 mm−1 Angiogenic vasculature development around a tumor site typically changes from the healthy vessel texture spacing of around 100 μm; k=10 mm−1 Due to its chaotic nature, the spacings in angiogenic vasculature cover a broad range from about 10 μm to 1 mm, or 1 mm−1 to 100 mm−1 Diagnostic assessment of dementia-related changes to the cortical neuronal spacing involves measuring high k-values, the healthy structure being about 100 μm spacing or k=10 mm−1 Variations of about 10-20% of this value, with increasing randomness in structure, mark the disease.
Therefore, the acquisition parameters can be chosen such that (1) the gradient height/duration generates a range of k-encodes spanning the neighborhood of k-space over which it is desired to inspect the power density present in the targeted tissue texture, (2) the samples to be combined must occur close enough in time that there is no significant blurring due to subject motion across the acquisition time of a block of samples to be combined, and (3) the tolerable amount of motion depends on the neighborhood of k-space under investigation (i.e., the wavelength).
Acquisition of textural data from within a targeted VOI with the non-zero gradient enables determination of the local variation of power density of k-values within a neighborhood of the initial k-value in k-space. The extent of k-space sampled at each gradient pulse is determined by the gradient height and the pulse width. Spacing between the samples in k-space is determined by the gradient height and the sampling rate.
These parameters are selected (1) to allow acquisition of sufficient data for combining toward significant SNR improvement, before subject motion can blur the data significantly relative to the texture to be measured, (2) to ensure sufficient correlation across the blocks of k-values from the acquisition to be combined to maintain a SNR≧0.5 dB, and (3) to set the extent of k-space over which the power density of k-values present in the texture is desired.
Blocks of sequential signal samples to be recombined for SNR improvement can be non-overlapping, or overlapping by a selected number of points, or a sliding block used so as to combine, for example, measures 1-4, 2-5, 3-6 and so on as will be described subsequently. Additionally, the number of samples in each block may be varied from block to block across the extent in k-space of the acquisition, this variation in number of samples to be combined being determined by the requirement for sufficient correlation to maintain SNR sufficient to provide a robust measurement. The approximate noise level can be determined independently by several methods well known to the industry including measuring noise in the absence of signal input.
Acquiring data with different magnitude gradients within one echo, TR, or scan may be accomplished with the successive gradient heights being selected to enable best SNR of the combined signal at the various targeted regions of k-space. To enable combination of sequential samples from the ADC to improve SNR, correlation among successive samples can be increased by proper choice of windowing in profile space, a shorter window driving a greater correlation distance across sequential values in k-space and a longer window resulting in lower sample to sample correlation. The window width selected is defined by both the desire for correlation across many samples in k-space, which dictates a shorter window, and the need to sample a sufficient extent of texture in real space to provide robust measure, especially when measuring highly amorphous textures.
Post-acquisition combination of the signals acquired in k-space in blocks, the number of samples to be combined determined by the requirement that the correlation between the individual signals to be combined be sufficient to achieve a SNR≧0 dB (the level of correlation is determined by subject motion, gradient height, sampling rate, window shape, and the underlying texture.)
Use of non-zero gradient acquisition may be employed to intentionally vary the direction and magnitude of the k-vector over a range during data acquisition, for the purpose of smoothing signal speckle—which will manifest as a time varying signal during the data acquisition—that results from interference of the varying phases and amplitudes of the individual spin signals. The selected variation in k-value direction and magnitude during data acquisition is chosen to provide sufficient combined measures to get a an estimation of the representative power within a neighborhood of k-space, with a SNR of 0 dB, where the neighborhood is within 20% of the 3D orientation and magnitude of the centroid of the neighborhood.
Correction for variation in signal magnitude at k-values across a set neighborhood created by application of a non-zero gradient may be accomplished by employing prescribed k-value encodes for a specific set of signal magnitude measurements at those k-values. Additionally, correlation within a set of measurements at k-values acquired within a time period and from a selected VOI can be induced by selecting the time period such that the biological motion is sufficiently small that the phase shift in the data induced by patient motion is less than 50% of the wavelength corresponding to the targeted textural k-value range. Alternatively, a windowing function may be selected such that there is sufficient correlation between individual measurements and the set estimate that a desired SNR can be achieved.
Details of the very-low SNR acquisition mode at even further expanded scale are shown in FIG. 13. In this portion of the sequence, the k-value is constant at an initial value 1114a, a second value 1114b induced by k value selection gradient pulse 1100a, a third value 1114c induced by k value selection gradient pulse 1100b and a fourth value 1114d induced by k value selection gradient pulse 1100c in the region 1116 where the sample gate is open thereby producing samples 1118. This is the previously described pulse sequence where multiple repeats of signal at the same k value are rapidly sampled, all of which are then combined into one estimate.
Details of the low SNR acquisition at the further expanded scale are shown in FIG. 14. In this portion of the sequence, notice that the k-values do change as shown by trace segment 1120, albeit slowly, due to the non-zero time-dependent phase-encode gradient 1106 present during the recording of the region 1122 when the sample gate is open. However, the range of samples 1124 across k-space is a relatively compact neighborhood where the values are highly correlated due to the windowing function.
Details of the high SNR acquisition at the further expanded scale are shown in FIG. 15. In this portion of the sequence, the k-values again change as shown by trace segment 1126, due to the non-zero time-dependent phase-encode gradient 1112 present during the opening of the sample gate in region 1128. The range of samples 1130 across K-space is still a relatively compact neighborhood, but outside the inter-sample correlation imposed by the window function.
The low SNR and high SNR acquisition modes with non-zero gradient are distinct from a standard frequency-encoded MRI sequence as the applied gradient is not used to establish a position, i.e. frequency encoding, but as a time dependent phase encode to rapidly acquire a number of individual samples across a relatively broader neighborhood of k-space.
As previously asserted, gradient acquisition can be acquired in one echo rather than in multiple echoes as seen FIG. 16. Again the illustrated pulse sequence for selecting the desired VOI and initial phase wrap for k is as described in FIG. 3 and is numbered consistently in FIG. 16.
As seen in FIG. 16 and at larger scale in FIG. 17, the k-value is constant at an initial value 1614a, a second value 1614b induced by k value selection gradient pulse 1600a and a third value 1614c induced by k value selection gradient pulse 1600b. Note that the k values are decremented as opposed to incremented in the example of FIG. 12. This is again the previously described pulse sequence where multiple repeats of signal at the same K value are rapidly sampled 1622a, 1622b and 1622c, all of which are then combined into one estimate.
In a second portion of the sequence, within the same echo, the k-values do change as shown by trace segment 1620, albeit slowly, due to the non-zero time-dependent phase-encode gradient 1606 present during sampling. However, the range of samples 1624 across k-space is a relatively compact neighborhood where the values are highly correlated.
In a third portion of the sequence, again still within the same echo, high SNR acquisition is conducted. The k-values again change as shown by trace segment 1626, due to the non-zero time-dependent phase-encode gradient 1608 present during the opening of the sample gate. The range of 1628 across k-space is still a relatively compact neighborhood but outside the inter-sample correlation imposed by the window function.
As seen in FIG. 18 (and in larger scale in FIG. 19), where the illustrated pulse sequence for selecting the desired VOI and initial phase wrap for k is as described in FIG. 3 and is numbered consistently in FIG. 18, a low non-zero magnitude gradient 1802 acting as a time dependent phase encode is applied and data samples 1804 are taken from an initial k-value 1806 for slowly time varying k-values, seen in trace segment 1808, having high correlation as previously described. The initial phase wrap may be selected to provide an initial k-value with a magnitude corresponding to a low SNR region.
Similarly, as seen in FIG. 20 (and in larger scale in FIG. 21), a pulse sequence for selecting the desired VOI and initial phase wrap to set the k-value region is illustrated and is as described in FIG. 3 and is numbered consistently in FIG. 20. A higher non-zero gradient 2002 acting as a time dependent phase encode is applied and data samples 2004 are taken from an initial k-value 2006 for more rapidly time varying k-values, seen in trace segment 2008. The initial phase wrap may be selected to provide an initial k-value with a magnitude corresponding to a higher SNR region. The encoding gradient 326 may be employed to wind up to the lowest or highest k-value in a targeted texture and the non-zero magnitude gradient pulse is imposed in the necessary direction (increasing or decreasing k) to reach the other limit in k-space to define the texture.
The acquired samples may be outside the inter-sample correlation imposed by the window function. However, signal levels for the k-values are relatively large and have high SNR. Additionally as previously described, rapid acquisition of samples in subsets 2010a, 2010b, 2010c and 2010d as exemplary, may be accomplished in a manner that the samples within the subset may remain sufficiently correlated and may provide desired data in structures having predetermined or anticipated texture. One can combine however many sequential values are correlated enough to yield an improvement in SNR via the combination (averaging being one simple form of combining). Then, the set of combined data points is used to characterize the power distribution across the entire acquisition, to get a better measure of the underlying texture within the VOI.
As previously discussed, rephrasing separate low and high k values based on low phase change in a second 90-180-180 excitation (TR). SNR is maximized with gradient ON acquisition by smart combination of successive k-value samples through reregistration of successive acquired signals. Data is acquired across a range of k-space for which the wavelengths are sufficiently long that subject motion can be easily corrected for by reregistration—i.e. the phase shift induced in the measure in this k range is much less than the textural wavelength. The low k-value signal is sampled in alternate refocusing sequences, or sequential excitations (TR), with the acquisition of the signal from the higher k-value range of interest. TE long wavelength measure is used to determine the motion-induced phase shift across the measurements. That phase shift is then applied to the higher k-data prior to reregistration.
A number of correlations are implied by spatial windowing.
If g(x) corresponds to a 1D (real-valued) signal, the corresponding function in K-space is given by the Fourier transform as:

G(2πk)=∫−∞∞g(x)e−j2πxkdx  (1)

Which is frequently expressed as a Fourier pair as

g(x)G(2πk)  (2)

Windowing is the process of limiting the extent of g(x) to a finite region of compact support, but doing it in such a way to minimize spectral artifacts due to discontinuities (artificially) introduced by the truncation.
Despite the specific shape used of the window function, there is an inverse relationship between the width of the window, and its spectrum. This is due to the Fourier relationship




  
    
      
        
          h
          ⁡
          
            (
            ax
            )
          
        
        ⇔
        
          
            1
            
              
              α
              
            
          
          ⁢
          
            H
            ⁡
            
              (
              
                
                  2
                  ⁢
                  π
                  ⁢
                  
                    
                  
                  ⁢
                  k
                
                α
              
              )
            
          
        
      
    
    
      
        (
        3
        )
      
    
  




Multiplying two functions has the effect of convolving their respective spectra, i.e.

f(x):=g(x)h(x)F(2πk)=G(2πk)*H(2πk)  (4)

The convolution can be thought of as a linear filtering of the spectrum as though the spectrum was the input signal
The term




  H
  ⁡
  
    (
    
      
        2
        ⁢
        π
        ⁢
        
          
        
        ⁢
        k
      
      α
    
    )
  




acts like a low-pass filter to the G(2πk) spectrum, which tends to smooth out the signal: the larger the value of a, the narrower the low-pass filter. This creates a significant correlation between adjacent values of F(2πk).

Estimators which observe noisy samples of a filtered input are well studied and can be applied to generate optimal estimates; Weiner filters, Kalman filters, etc.
Dynamic acquisition modes may be employed wherein:
X corresponds to a 3D vector in image space,
g(X) corresponds to the value of the image at a given 3D spatial location,
K corresponds to the 3D vector,
G(K) corresponds to the value in k-space of the image g.
For initial simplicity, the time-dependency of this signal is ignored which in turn depends upon T1, T2, T2*, as well as signal contribution due to differing isochromats (different chemical species within the Volume) etc. In the sequel the effect of these is taken into account
Basic Principles relied upon are:
Generally, SNR of G(K) is highest at |k|=0, then decreases with increasing |k|
The rate at which SNR decreases is typically expressed as SNR ∝|k|−α where α is in the range of 1-3.
The sampling rate, combined with the magnitude of the gradient will set the sample spacing (Δk) density in k-space for a given VOI.
As the gradient magnitude is decreased, the sample density increases (i.e. Δk decreases). Depending upon the size of the windowing in image space, there is a corresponding correlation implied.
For a generalized case the simplified MRI relationship between spatial coordinates and K-space given by




  
    
      
        
          S
          ⁡
          
            (
            k
            )
          
        
        =
        
          
            
              ∫
              
                ∫
                ∫
              
            
            
              
              3
            
          
          ⁢
          
            I
            ⁡
            
              (
              r
              )
            
          
          ⁢
          
            e
            
              
                -
                j
              
              ⁢
              
                
              
              ⁢
              2
              ⁢
              π
              ⁢
              
                
              
              ⁢
              k
              ⁢
              
                
              
              ⁢
              •
              ⁢
              
                
              
              ⁢
              r
            
          
          ⁢
          
            d
            r
          
        
      
    
    
      
        (
        6
        )
      
    
  




Where

r represents the real valued 3-Dimensional spatial coordinates with units of meters (m).

I(r) represents the image which is a non-negative Real function of spatial coordinates r.

k represents the real valued 3-Dimensional k-space coordinates with units in cycles/meter (m−1)
S(k) represents the Fourier Transform of I(r) and is generally a complex-valued function of k
And the integral is over the entire 3-Dimensional spatial plane.
In words, S(k) represents the corresponding value in 3-dimensional k-space of the image function/(r).
The k-space coordinates, in turn, are a function of time and have the general form




  
    
      
        
          k
          ⁡
          
            (
            t
            )
          
        
        =
        
          
            γ
            _
          
          ⁢
          
            
              ∫
              
                -
                ∞
              
              t
            
            ⁢
            
              
                g
                ⁡
                
                  (
                  τ
                  )
                
              
              ⁢
              
                
              
              ⁢
              
                d
                τ
              
            
          
        
      
    
    
      
        (
        7
        )
      
    
  




Where

 is the proton gyromagnetic ratio with value 42.576 MHz/T

g(t) is a real-valued 3-Dimensional function of time representing the gradient strength with units in T/m. This function, is a design input as part of the pulse sequence whose purpose is to manipulate the proton spins in some desired way.

The integral in equation (7) indicates that the value of k(t) for a given value oft, is computed as the integral of all previous history of the gradient function. While technically correct, it is often more convenient to express this as




  
    
      
        
          k
          ⁡
          
            (
            t
            )
          
        
        =
        
          
            
              γ
              _
            
            ⁢
            
              
                ∫
                
                  t
                  0
                
                t
              
              ⁢
              
                
                  g
                  ⁡
                  
                    (
                    τ
                    )
                  
                
                ⁢
                
                  
                
                ⁢
                
                  d
                  τ
                
              
            
          
          +
          
            (
            
              t
              0
            
            )
          
        
      
    
    
      
        (
        8
        )
      
    
  




Where now t0 represents a convenient starting time, k(t0) is the corresponding k-value at t0, and the lower limit of the integral starts at t0.
Making the dependence on time more explicit, equation (6) can be expressed as




  
    
      
        
          S
          ⁡
          
            (
            t
            )
          
        
        =
        
          
            
              
                ∫
                
                  ∫
                  ∫
                
              
              
                
                3
              
            
            ⁢
            
              I
              ⁡
              
                (
                r
                )
              
            
            ⁢
            
              e
              
                
                  -
                  j
                
                ⁢
                
                  
                
                ⁢
                2
                ⁢
                
                  
                
                ⁢
                π
                ⁢
                
                  
                
                ⁢
                
                  k
                  ⁡
                  
                    (
                    t
                    )
                  
                
                ⁢
                r
              
            
            ⁢
            
              d
              r
            
            ⁢
            
              
            
            ⁢
            t
          
          ≥
          
            t
            0
          
        
      
    
    
      
        (
        9
        )
      
    
  




S(t) represents the complex-valued baseband signal one might obtain during an MRI echo experiment which is played in conjunction with a gradient sequence encoded ing(t).

Without loss of generality, k, g and r can be decomposed into Cartesian components as

k=[kxkykz]T 

g=[gxgygz]T 

r=[rxryrz]T  (10)

And express (9) as





  
    
      
        
          
            S
            ⁡
            
              (
              t
              )
            
          
          =
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                ∫
                
                  -
                  ∞
                
                ∞
              
              ⁢
              
                
                  ∫
                  
                    -
                    ∞
                  
                  ∞
                
                ⁢
                
                  
                    I
                    ⁡
                    
                      (
                      
                        
                          r
                          x
                        
                        ,
                        
                          r
                          y
                        
                        ,
                        
                          r
                          z
                        
                      
                      )
                    
                  
                  ⁢
                  
                    e
                    
                      
                        -
                        j
                      
                      ⁢
                      
                        
                      
                      ⁢
                      2
                      ⁢
                      
                        
                      
                      ⁢
                      
                        π
                        ⁡
                        
                          (
                          
                            
                              
                                
                                  k
                                  x
                                
                                ⁡
                                
                                  (
                                  t
                                  )
                                
                              
                              ⁢
                              
                                r
                                x
                              
                            
                            +
                            
                              
                                
                                  k
                                  y
                                
                                ⁡
                                
                                  (
                                  t
                                  )
                                
                              
                              ⁢
                              
                                r
                                y
                              
                            
                            +
                            
                              
                                
                                  k
                                  z
                                
                                ⁡
                                
                                  (
                                  t
                                  )
                                
                              
                              ⁢
                              
                                r
                                z
                              
                            
                          
                          )
                        
                      
                    
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    d
                    
                      r
                      x
                    
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    d
                    
                      r
                      y
                    
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    d
                    
                      r
                      z
                    
                  
                
              
            
          
        
        ⁢
        
          

        
        ⁢
        
          
        
        ⁢
        
          t
          ≥
          
            t
            0
          
        
      
    
    
      
        (
        11
        )
      
    
  




In general k(t) represents a curvilinear path within K-space as a function of time.
Initially, to facilitate explaining the initial concept, evaluation is confined along a single dimension by assuming k(t)=[kx(t) 0 0]T. Equation (11) then simplifies to




  
    
      
        
          
            S
            ⁡
            
              (
              t
              )
            
          
          =
          
            
              
                ∫
                
                  -
                  ∞
                
                ∞
              
              ⁢
              
                
                  ρ
                  ⁡
                  
                    (
                    
                      r
                      x
                    
                    )
                  
                
                ⁢
                
                  e
                  
                    
                      -
                      j
                    
                    ⁢
                    
                      
                    
                    ⁢
                    2
                    ⁢
                    
                      
                    
                    ⁢
                    π
                    ⁢
                    
                      
                    
                    ⁢
                    
                      
                        k
                        x
                      
                      ⁡
                      
                        (
                        t
                        )
                      
                    
                    ⁢
                    
                      r
                      x
                    
                  
                
                ⁢
                
                  
                
                ⁢
                
                  d
                  
                    r
                    x
                  
                
                ⁢
                
                  
                
                ⁢
                t
              
            
            ≥
            
              t
              0
            
          
        
        ⁢
        
          

        
        ⁢
        Where
      
    
    
      
        (
        12
        )
      
    
  
  
    
      
        
          ρ
          ⁡
          
            (
            
              r
              x
            
            )
          
        
        =
        
          
            ∫
            
              -
              ∞
            
            ∞
          
          ⁢
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                I
                ⁡
                
                  (
                  
                    
                      r
                      x
                    
                    ,
                    
                      r
                      y
                    
                    ,
                    
                      r
                      z
                    
                  
                  )
                
              
              ⁢
              
                
              
              ⁢
              
                d
                
                  r
                  y
                
              
              ⁢
              
                
              
              ⁢
              
                d
                
                  r
                  z
                
              
            
          
        
      
    
    
      
        (
        13
        )
      
    
  




And equation (8) reduces to





  
    
      
        
          
            k
            x
          
          ⁡
          
            (
            t
            )
          
        
        =
        
          
            
              γ
              _
            
            ⁢
            
              
                ∫
                
                  t
                  0
                
                t
              
              ⁢
              
                
                  
                    g
                    x
                  
                  ⁡
                  
                    (
                    τ
                    )
                  
                
                ⁢
                
                  
                
                ⁢
                
                  d
                  τ
                
              
            
          
          +
          
            
              k
              x
            
            ⁡
            
              (
              
                t
                0
              
              )
            
          
        
      
    
    
      
        (
        14
        )
      
    
  




Define





  
    
      
        
          R
          ⁡
          
            (
            k
            )
          
        
        :=
        
          
                        ⁢
            
              {
              
                ρ
                ⁡
                
                  (
                  x
                  )
                
              
              }
            
            ⁢
            
              (
              k
              )
            
          
          =
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                ρ
                ⁡
                
                  (
                  x
                  )
                
              
              ⁢
              
                e
                
                  
                    -
                    j
                  
                  ⁢
                  
                    
                  
                  ⁢
                  2
                  ⁢
                  
                    
                  
                  ⁢
                  π
                  ⁢
                  
                    
                  
                  ⁢
                  kx
                
              
              ⁢
              
                
              
              ⁢
              
                d
                x
              
            
          
        
      
    
    
      
        (
        15
        )
      
    
  




Which is occasionally expressed as

R(k)ρ(x)  (16)

To concisely indicate that R(k) and ρ(x) are Fourier Transform pairs. By comparing (12) to (15), it can be seen that S(t) is just a time dependent progression across various Fourier coefficients represented by

S(t)=R(k(t))  (17)

Where the mapping between the time-value t and the corresponding K-space coordinate is given by equation (14).

To generally model the receive signal In an actual MRI machine, a combination of the desired signal and noise received by the antenna. That signal is then filtered, amplified, down-converted, sampled, and quantized.
The specific details are machine-dependent, but a simple model can be developed to represent the output of the machine as follows:
Let Y(t) represent a combination of the signal of interest, and a noise signal as

Y(t)=s(t)+w(t)  (18)

Where

S(t)i is given in equation (17) and

W(t) is a complex-valued zero-mean, Additive White Gaussian Noise Process with variance σw2, i.e. E{W(t)}=0 and E{W(t)W*(t+τ)}=σw2δ(τ)

The received signal Y(t) is then uniformly sampled

Yn=Y(t)|t=n·Δt=R(k(t))|t=n·Δt+W(t)|t=n·Δt  (19)

Which can be expressed more simply as

Yn=R(kn)+Wn  (20)

Where the sequence kn is given by





  
    
      
        
          
            
              
                k
                0
              
              =
              
                k
                ⁡
                
                  (
                  t
                  )
                
              
            
            
          
          
            t
            =
            0
          
        
        ⁢
        
          

        
        ⁢
        
          
            k
            
              n
              +
              1
            
          
          =
          
            
              k
              n
            
            +
            
              
                γ
                _
              
              ⁢
              
                
                  ∫
                  
                    n
                    ⁢
                    
                      
                    
                    ⁢
                    Δ
                    ⁢
                    
                      
                    
                    ⁢
                    t
                  
                  
                    
                      (
                      
                        n
                        +
                        1
                      
                      )
                    
                    ⁢
                    Δ
                    ⁢
                    
                      
                    
                    ⁢
                    t
                  
                
                ⁢
                
                  
                    g
                    ⁡
                    
                      (
                      τ
                      )
                    
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    d
                    τ
                  
                
              
            
          
        
      
    
    
      
        (
        21
        )
      
    
  




Define




  
    
      
        
          Δ
          ⁢
          
            
          
          ⁢
          
            k
            
              n
              +
              1
            
          
        
        =
        
          
            γ
            _
          
          ⁢
          
            
              ∫
              
                n
                ⁢
                
                  
                
                ⁢
                Δ
                ⁢
                
                  
                
                ⁢
                t
              
              
                
                  (
                  
                    n
                    +
                    1
                  
                  )
                
                ⁢
                Δ
                ⁢
                
                  
                
                ⁢
                t
              
            
            ⁢
            
              
                g
                ⁡
                
                  (
                  τ
                  )
                
              
              ⁢
              
                
              
              ⁢
              
                d
                τ
              
            
          
        
      
    
    
      
        (
        22
        )
      
    
  




Then (21) can be simply expressed as

k0=k(t)|t=0 

kn+1=kn+Δkn+1  (23)

In words, then, the sequence kn is defined by a sequence of increments which is determined by the integral between samples of the gradient function.
Equations (20), (22) and (23) may be employed to describe the signals under different gradient conditions disclosed herein.
Collecting samples of an echo which has been “pre-phased” through some gradient activity before-hand, but now the gradient is no longer held to zero as described above with respect to FIGS. 12, 14 and 15 can be analyzed as follows.
The signal is then given by equation (20) as

Yn=R(kn)+Wn  (24)

And kn is given by equation (21) as





  
    
      
        
          
            
              
                k
                0
              
              =
              
                k
                ⁡
                
                  (
                  t
                  )
                
              
            
            
          
          
            t
            =
            0
          
        
        ⁢
        
          

        
        ⁢
        
          
            k
            
              n
              +
              1
            
          
          =
          
            
              k
              n
            
            +
            
              
                γ
                _
              
              ⁢
              
                
                  ∫
                  
                    n
                    ⁢
                    
                      
                    
                    ⁢
                    Δ
                    ⁢
                    
                      
                    
                    ⁢
                    t
                  
                  
                    
                      (
                      
                        n
                        +
                        1
                      
                      )
                    
                    ⁢
                    Δ
                    ⁢
                    
                      
                    
                    ⁢
                    t
                  
                
                ⁢
                
                  
                    g
                    ⁡
                    
                      (
                      τ
                      )
                    
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    d
                    τ
                  
                
              
            
          
        
      
    
    
      
        (
        25
        )
      
    
  




Since, measurement is occurring in a non-zero gradient regime, the integral term is no longer zero, which implies that the sequence kn is no longer constant, and in turn the sequence R(kn) is no longer constant.
Since no assumptions have been made on the underlying structure of I(r), it cannot be implied that there is any particular structure or relationship amongst the values of R(kn). This puts us at a distinct disadvantage when wanting to estimate useful signals in a very low SNR environment.
A structure may be imposed on the values of R(kn) by applying a multiplicative window function in the image domain. This is accomplished by leveraging two Fourier Transform identities:
Multiplication in one domain corresponds to convolution in the reciprocal domain.
Define the following Fourier Pairs:

ν(x)N(k)

ρ(x)R(k)

ζ(x)Z(k)  (26)

Then, the product in one domain corresponds to convolution in the reciprocal domain:

ν(x)=ρ(x)ζ(x)N(k)=R(k)*Z(k)  (27)

Scaling in one domain corresponds to an inverse scaling in the reciprocal domain.
If ζ(x)Z(k) then




  
    
      
        
          ζ
          ⁡
          
            (
            
              x
              a
            
            )
          
        
        ⇔
        
          
            
            a
            
          
          ⁢
          
            Z
            ⁡
            
              (
              
                a
                ·
                k
              
              )
            
          
        
      
    
    
      
        (
        28
        )
      
    
  




Windowing functions are typically used to limit the image space to a finite, compact region of interest, while at the same time, minimizing the adverse consequences on the corresponding image spectrum due to the window itself. Those skilled in the art will appreciate there are a wide variety of window functions which have been developed, each of which have their own particular set of characteristics.
For sake of illustration, consider the most basic window function:




  
    
      
        
          rect
          ⁡
          
            (
            t
            )
          
        
        =
        
          
            Π
            ⁡
            
              (
              t
              )
            
          
          =
          
            {
            
              
                
                  0
                
                
                  
                    
                      
                      t
                      
                    
                    >
                    
                      1
                      2
                    
                  
                
              
              
                
                  
                    1
                    2
                  
                
                
                  
                    
                      
                      t
                      
                    
                    =
                    
                      1
                      2
                    
                  
                
              
              
                
                  1
                
                
                  
                    
                      
                      t
                      
                    
                    <
                    
                      1
                      2
                    
                  
                
              
            
          
        
      
    
    
      
        (
        29
        )
      
    
  




The corresponding Fourier transform is given by





  
    
      
        
          F
          ⁢
          
            {
            
              Π
              ⁡
              
                (
                t
                )
              
            
            }
          
        
        =
        
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                Π
                ⁡
                
                  (
                  t
                  )
                
              
              ⁢
              
                e
                
                  
                    -
                    j
                  
                  ⁢
                  
                    
                  
                  ⁢
                  2
                  ⁢
                  
                    
                  
                  ⁢
                  π
                  ⁢
                  
                    
                  
                  ⁢
                  ft
                
              
              ⁢
              
                
              
              ⁢
              
                d
                t
              
            
          
          =
          
            
              
                sin
                ⁡
                
                  (
                  
                    π
                    ⁢
                    
                      
                    
                    ⁢
                    f
                  
                  )
                
              
              
                π
                ⁢
                
                  
                
                ⁢
                f
              
            
            =
            
              sin
              ⁢
              
                
              
              ⁢
              
                c
                ⁡
                
                  (
                  f
                  )
                
              
            
          
        
      
    
    
      
        (
        30
        )
      
    
  




which is frequently expressed as the Fourier pair

Π(t)sin c(f)  (31)

Using equation (28) a slightly generalized version and its Fourier pair is





  
    
      
        
          Π
          ⁡
          
            (
            
              t
              T
            
            )
          
        
        ⇔
        
          
            
            T
            
          
          ⁢
          sin
          ⁢
          
            
          
          ⁢
          
            c
            ⁡
            
              (
              Tf
              )
            
          
        
      
    
    
      
        (
        32
        )
      
    
  




Using equation (27), the windowed profile and Fourier pair is





  
    
      
        
          v
          ⁡
          
            (
            x
            )
          
        
        =
        
          
            
              
                Π
                ⁡
                
                  (
                  
                    x
                    X
                  
                  )
                
              
              ⁢
              
                ρ
                ⁡
                
                  (
                  x
                  )
                
              
            
            ⇔
            
              N
              ⁡
              
                (
                k
                )
              
            
          
          =
          
            
              
              X
              
            
            ⁢
            sin
            ⁢
            
              
            
            ⁢
            
              c
              ⁡
              
                (
                Xk
                )
              
            
            *
            
              R
              ⁡
              
                (
                k
                )
              
            
          
        
      
    
    
      
        (
        33
        )
      
    
  




Using equation (24) as a reference, the sampled MRI signal can be expressed as

Yn=N(kn)+Wn  (34)

Which, using (33), can be expanded as





  
    
      
        
          Y
          n
        
        =
        
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                
                X
                
              
              ⁢
              sin
              ⁢
              
                
              
              ⁢
              
                c
                ⁡
                
                  (
                  Xq
                  )
                
              
              ⁢
              
                R
                ⁡
                
                  (
                  
                    
                      k
                      n
                    
                    -
                    q
                  
                  )
                
              
              ⁢
              
                
              
              ⁢
              
                d
                q
              
            
          
          +
          
            W
            n
          
        
      
    
    
      
        (
        35
        )
      
    
  




Where the convolution integral has been specifically expanded.

The value of the convolution integral taken at kn is no longer a function of just one point of R(kn). For each point kn the convolution integral computes a weighted sum of the values of R(k) centered around kn The extent of the neighborhood in k-space is inversely proportional to the parameter X: Smaller values of X increase the width of the neighborhood in k-space.
For the embodiments herein the extent of the domain of values of interest correspond to the collection of k-space values k0, k1, k2, . . . kN-1. Define




  
    
      
        
          
            k
            min
          
          =
          
            
              min
              n
            
            ⁢
            
              k
              n
            
          
        
        ⁢
        
          

        
        ⁢
        
          
            k
            max
          
          =
          
            
              max
              n
            
            ⁢
            
              k
              n
            
          
        
      
    
    
      
        (
        36
        )
      
    
  




Which in turn are functions of the time interval Δt and the function g(τ).

For example, making a simplifying assumption that g(τ)=G where G is a positive constant, then kn is just a uniform sampling across a portion of k-space, and is given by

kn=k0+nGΔt  (37)

Then kmin and kmax is given by

kmin=k0 

kmax=k0+(N−1)GΔt  (38)

While a simple sampling of k-space may be chosen, it is not specifically required. Indeed there could be applications where non-uniform and/or even non-monotonic sampling strategies could be useful.
Ideally, the parameter X (and the window function) are chosen so that the resulting weighted sum across the neighborhood of kn is “wide enough” so that N(kn)≈C
where C is a complex-valued constant, but not so wide so as to lose significant spectral resolution
For purposes of the disclosed embodiments herein a “small” non-zero gradient may be determined based on selection of desired windowing. From equation (10)




  
    
      
        
          R
          ⁡
          
            (
            k
            )
          
        
        :=
        
          
                        ⁢
            
              {
              
                ρ
                ⁡
                
                  (
                  x
                  )
                
              
              }
            
            ⁢
            
              (
              k
              )
            
          
          =
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                ρ
                ⁡
                
                  (
                  x
                  )
                
              
              ⁢
              
                e
                
                  
                    -
                    j
                  
                  ⁢
                  
                    
                  
                  ⁢
                  2
                  ⁢
                  
                    
                  
                  ⁢
                  π
                  ⁢
                  
                    
                  
                  ⁢
                  kx
                
              
              ⁢
              
                
              
              ⁢
              
                d
                x
              
            
          
        
      
    
    
      
        (
        41
        )
      
    
  




Assume that the nominal center point of the profile has shifted to be centered around a point x0. This results in




  
    
      
        
          
            R
            
              x
              0
            
          
          ⁡
          
            (
            k
            )
          
        
        :=
        
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                ρ
                ⁡
                
                  (
                  
                    x
                    -
                    
                      x
                      0
                    
                  
                  )
                
              
              ⁢
              
                e
                
                  
                    -
                    j
                  
                  ⁢
                  
                    
                  
                  ⁢
                  2
                  ⁢
                  
                    
                  
                  ⁢
                  π
                  ⁢
                  
                    
                  
                  ⁢
                  kx
                
              
              ⁢
              
                
              
              ⁢
              
                d
                x
              
            
          
          =
          
            
              e
              
                
                  -
                  j
                
                ⁢
                
                  
                
                ⁢
                2
                ⁢
                
                  
                
                ⁢
                π
                ⁢
                
                  
                
                ⁢
                
                  kx
                  0
                
              
            
            ⁢
            
              R
              ⁡
              
                (
                k
                )
              
            
          
        
      
    
    
      
        (
        42
        )
      
    
  




Which indicates that each point in k-space is rotated in complex space proportional to the offset x0.

It can be assumed that the gradient is a positive constant, then, by equation (39)

kn=k0+nGΔt  (43)

Substituting in (42) produces

Rx0(kn)=e−jθ0e−jnΔθR(kn)  (44)

Where the initial phase offset θ0 and the phase increment Δθ is given by

θ0:=−2πk0x0 

Δθ:=−2πGΔtx0  (45)

In the event that, due to the application of a properly specified windowing function, R(kn)≈C a complex constant within the neighborhood, the post-acquisition estimator would first multiply an offsetting phase increment ejnΔθ to each acquired sample of Rx0(kn) before combining and generating the final estimate.
An estimate of Δθ can be obtained from a sequence of k-space samples taken of the windowed profile over lower k-values (where the SNR is higher).
Correlation may be induced by windowing as one parameter as discussed previously. Multiplication of the profile ρ(x) by a real-valued window function ζ(x) corresponds to convolution in k-space by the Fourier relation

ν(x)=ρ(x)ζ(x)N(k)=R(k)*Z(k)  (46)

Z(k) is treated as an impulse response of a linear filter which is applied to the complex-valued signal R(k) in k-space to produce a complex-valued output signal N(k).
The autocorrelation function of the output signal RNN(κ1,κ2) can be expressed as a function the autocorrelation of the input signal RRR(κ1,κ2), and the impulse response Z(k) as




  
    
      
        
          
            R
            NN
          
          ⁡
          
            (
            
              
                κ
                1
              
              ,
              
                κ
                2
              
            
            )
          
        
        =
        
          
            ∫
            
              -
              ∞
            
            ∞
          
          ⁢
          
            
              ∫
              
                -
                ∞
              
              ∞
            
            ⁢
            
              
                
                  R
                  RR
                
                ⁡
                
                  (
                  
                    
                      
                        κ
                        1
                      
                      -
                      α
                    
                    ,
                    
                      
                        κ
                        2
                      
                      -
                      β
                    
                  
                  )
                
              
              ⁢
              
                Z
                ⁡
                
                  (
                  α
                  )
                
              
              ⁢
              
                Z
                ⁡
                
                  (
                  β
                  )
                
              
              ⁢
              
                d
                α
              
              ⁢
              
                
              
              ⁢
              
                d
                β
              
            
          
        
      
    
    
      
        (
        47
        )
      
    
  




Equation (47) is inconvenient because the autocorrelation function of the underlying signal RRR(κ1,κ2) is not usually known. A simplifying assumption is made that R(k) is a white-noise, wide-sense stationary process, and express the autocorrelation as

RRR(κ1,κ2)=σR2RZZ(κ1−κ2)  (48)

With this assumption, (47) reduces to

RNN(κ1,κ2)×σR2RZZ(κ2−κ1)  (49)

Where RZZ(κ) is the autocorrelation function of the impulse response Z(k) and is given by





  
    
      
        
          
            R
            ZZ
          
          ⁡
          
            (
            κ
            )
          
        
        :=
        
          
            ∫
            
              -
              ∞
            
            ∞
          
          ⁢
          
            
              Z
              ⁡
              
                (
                k
                )
              
            
            ⁢
            
              Z
              ⁡
              
                (
                
                  k
                  +
                  κ
                
                )
              
            
            ⁢
            
              d
              k
            
          
        
      
    
    
      
        (
        50
        )
      
    
  




The mapping of k vs time is mapped as follows.

kn=k0+nGΔt  (51)

The normalized correlation





  
    
      
        
          η
          ⁡
          
            (
            
              NG
              ⁢
              
                
              
              ⁢
              Δ
              ⁢
              
                
              
              ⁢
              t
            
            )
          
        
        :=
        
          
            
              
                R
                NN
              
              ⁡
              
                (
                
                  
                    k
                    0
                  
                  ,
                  
                    
                      k
                      0
                    
                    +
                    
                      NG
                      ⁢
                      
                        
                      
                      ⁢
                      Δ
                      ⁢
                      
                        
                      
                      ⁢
                      t
                    
                  
                
                )
              
            
            
              
                R
                NN
              
              ⁡
              
                (
                
                  
                    k
                    0
                  
                  ,
                  
                    k
                    0
                  
                
                )
              
            
          
          =
          
            
              
                R
                ZZ
              
              ⁡
              
                (
                
                  NG
                  ⁢
                  
                    
                  
                  ⁢
                  Δ
                  ⁢
                  
                    
                  
                  ⁢
                  t
                
                )
              
            
            
              
                R
                ZZ
              
              ⁡
              
                (
                0
                )
              
            
          
        
      
    
    
      
        (
        4
        )
      
    
  




measures the degree to which the underlying sample points are correlated. In low SNR regimes, a high correlation is desired across all of the samples and therefore establish a lower bound:





  
    
      
        
          η
          min
        
        ≤
        
          
            
              R
              ZZ
            
            ⁡
            
              (
              
                NG
                ⁢
                
                  
                
                ⁢
                Δ
                ⁢
                
                  
                
                ⁢
                t
              
              )
            
          
          
            
              R
              ZZ
            
            ⁡
            
              (
              0
              )
            
          
        
      
    
    
      
        (
        53
        )
      
    
  




Equations (50) and (53) provide the defining relationship between the window function impulse response Z(k), the gradient strength G, the sample interval Δt, the number of samples N, and the correlation lower bound ηmin.
For example, assume that the window function is defined to be




  
    
      
        
          ζ
          ⁡
          
            (
            x
            )
          
        
        =
        
          Π
          ⁡
          
            (
            
              x
              X
            
            )
          
        
      
    
    
      
        (
        54
        )
      
    
  




Where Π(x) is a standard so-called rectangular function defined below, and X is a constant.





  
    
      
        
          Π
          ⁡
          
            (
            x
            )
          
        
        =
        
          {
          
            
              
                0
              
              
                
                  
                    
                    x
                    
                  
                  >
                  
                    1
                    2
                  
                
              
            
            
              
                
                  1
                  2
                
              
              
                
                  
                    
                    x
                    
                  
                  =
                  
                    1
                    2
                  
                
              
            
            
              
                1
              
              
                
                  
                    
                    x
                    
                  
                  <
                  
                    1
                    2
                  
                
              
            
          
        
      
    
    
      
        (
        55
        )
      
    
  




The impulse function Z(k) is given by the Fourier transform




  
    
      
        
          Z
          ⁡
          
            (
            k
            )
          
        
        =
        
          
            
              
              X
              
            
            ⁢
            
              
                sin
                ⁡
                
                  (
                  
                    π
                    ⁢
                    
                      
                    
                    ⁢
                    Xk
                  
                  )
                
              
              
                π
                ⁢
                
                  
                
                ⁢
                Xk
              
            
          
          =
          
            
              
              X
              
            
            ⁢
            sin
            ⁢
            
              
            
            ⁢
            
              c
              ⁡
              
                (
                Xk
                )
              
            
          
        
      
    
    
      
        (
        56
        )
      
    
  




The corresponding normalized correlation function η(κ) is given by

η(κ)=sin c(Xκ)  (57)

Restricting the correlation to be lower bounded by ηmin=0.95 then, by (53) the condition arises that

ηmin≦sin c(X·N·G·Δt)  (58)

This can be approximated using the first two terms of a Taylor series as





  
    
      
        
          
            sin
            ⁡
            
              (
              θ
              )
            
          
          θ
        
        ≈
        
          1
          -
          
            
              θ
              2
            
            
              3
              !
            
          
        
      
    
    
      
        (
        59
        )
      
    
  




Which can be inverted and applied to (58) to produce





  
    
      
        
          
            N
            ·
            G
            ·
            Δ
          
          ⁢
          
            
          
          ⁢
          t
        
        ≤
        
          
            
              
                3
                !
              
              ⁢
              
                (
                
                  1
                  -
                  
                    η
                    min
                  
                
                )
              
            
          
          
            π
            ⁢
            
              
            
            ⁢
            X
          
        
      
    
    
      
        (
        60
        )
      
    
  




Which now explicitly expresses an upper bound on the product of the gradient strength G, the sample interval Δt, and the number of samples N.

Typically the sample interval Δt, and the number of samples N are determined by other considerations. Taking these as given, the maximum gradient level is then given by




  
    
      
        G
        ≤
        
          
            
              
                3
                !
              
              ⁢
              
                (
                
                  1
                  -
                  
                    η
                    min
                  
                
                )
              
            
          
          
            π
            ⁢
            
              
            
            ⁢
            
              X
              ·
              N
              ·
              Δ
            
            ⁢
            
              
            
            ⁢
            t
          
        
      
    
    
      
        (
        61
        )
      
    
  




For a non-zero gradient data acquisition in this case, as long as the gradient G is below the calculated upper bound, the samples acquired will have the defined correlation level. This condition as defined for purposes herein as a “small gradient” level.
Sampling past this limitation will result in lower sample correlation and therefore have less of a potential post acquisition SNR gain. A “higher” gradient may be defined as operating in this condition. Gradient determination is affected by a number of parameters including (1) choice of the window function (e.g. rectangular, Tukey, Hamming, etc.) which influences the shape (and to a certain extent, the width) of the “main lobe” in the impulse response, (2) choice of window extent (the larger the extent in the profile domain, the narrower the “main lobe” in the impulse response), (3) the impulse response which may create an autocorrelation function, (4) the desired level of correlation which determines the effective width in k-space, within which the samples must be contained, and (5) sampling rate*Number of samples*gradient size which determines the actual sampling neighborhood size (note, as long as this number is bounded by the number contained in element (4)) the gradient remains in the “lower gradient level” regime.
An exemplary embodiment maintains a constant ratio of textural wavelength to length of VOI acquisition axis. As the targeted k-value varies, the length of the VOI acquisition axis is varied such that the ratio of the corresponding textural wavelength to the acquisition length remains constant. The aim here is to keep the number of textural “cells” sampled constant. In this way the differential broadening observed at specific points in k-space, Δk, is expected to arise from sources other than sampled length in real space, such as the finite width of the RF pulse or the edges of the gradient pulse.
MR-based diagnostic Techniques may be combined. Certain MR-based techniques designed to look at very fine tissue structure provide data that may be difficult to interpret in certain pathologies, as they provide only an indirect measure of the underlying structures. Diffusion weighted imaging and Magnetic Resonance Elastography (MRE) are two such techniques. The method of this provisional filing is a direct measure and hence would provide, in many cases, a better measure of fine texture, and in some cases provides complementary data to increase diagnostic capability. Combining acquisition techniques can provide more robust measure of texture, and hence of pathology.
The embodiments disclosed may be used in combination with Magnetic Resonance Elastography (MRE). Currently, the main application of MRE is as a diagnostic for liver disease to determine therapy response, progression, need for biopsy, etc. Though the targeted pathology is fibrotic development, the technique measures this indirectly, through measurement of tissue stiffness. In many cases, it is difficult to distinguish fibrotic development from other stiffness-inducing conditions such as portal hypertension and inflammation. Further, hepatic iron overload, which often results from a compromised liver, will lead to low signal, hence inadequate visualization of the induced mechanical waves.
The method disclosed herein can provide direct measure of fibrotic development in the liver and, as such, would provide additional data on disease progression or response to therapy in the case of the various triggers of fibrotic liver disease. It provides a local measure within the targeted anatomy for calibration of other, indirect measures, such as MRE, DTI, DWI, etc.
The embodiments disclosed may be used in combination with, or replacement for, diffusion weighted imaging in tumors. The ability to detect the edge of tumors with high accuracy would facilitate accurate surgical removal. Can acquire data using the method disclosed herein in VOIs along a direction through a tumor region looking for the edge of the region of angiogenic vasculature.
The ability to measure inside of tumors to gauge therapy response would help in targeting intervention. As an example of the latter, immunotherapy treatment of melanoma tumors induces swelling of the tumor due to T-cell infiltration which, on a structural MR scan, looks similar to malignant tumor growth. Hence, it is difficult to decide whether to continue the therapy. The ability to look at the state of the vasculature within the tumor would enable discernment of whether the growth was cancerous or due to immune system response. The method disclosed herein would also provide local calibration of the currently used DTI measures, which are often difficult to interpret.
The embodiments disclosed may be used as a bone degradation measure in oncology. It is well known that radiation and/or chemotherapy often compromise bone health. A measure of changes to bone resulting from cancer therapy would help in tailoring therapy and determine if there is need for intervention to protect bone health.
Currently, as a follow on to surgery and treatment for breast cancer, patients are routinely put in the MR scanner to image the breast tissue. The sternum is within the field of view for such exams, enabling easy application of a short add-on sequence of this method to measure changes to trabecular bone and thus obtain a measure of bone health.
As further examples of potential use of the method disclosed herein in oncology, the embodiments disclosed may be used to measure and quantify hyperplasiac development of mammary duct growth in response to tumor formation and development or to measure and quantify angiogenic growth of vasculature surrounding tumors to stage development, type, and response to therapy. Ongoing treatment after breast surgery often involves reducing estrogen levels, further compromising bone health and, as such, referral for MR scans for bone monitoring is common; use of the method of this patent filing would enable robust and detailed evaluation of bone health by direct measurement of the trabecular bone structure.
The disclosed embodiments are also complementary with Big Data and machine learning schemes. The method disclosed also complements the trend towards use of comparison among large aggregates of medical data to learn more about disease, increase predictive power for individual patients and for specific diseases, and note trends across various populations. Benefits of using the method of this filing in conjunction with Big data/machine learning include:
instead of 20 instances of the method of this patent output compared over a population of unknown pathology, for example, k-value distribution vs. fracture in femur may be evaluated, or changing k-values in cortical neuron bundles can be measured and correlated with performance on Alzheimer's mini-mental state exams or other assessment of AD, or the local k-values in liver compared with other inferences of liver disease over a huge population;

use of machine learning over large populations enables determination of specific biomarkers in pathology;

the ability to make useful correlations in big data gets much better with high SNR measure input such as that provided by the method of the embodiments disclosed;

such machine learning can indicate, for example, if a disease is defined by appearance of a specific k-value appearing in the diseased tissue.

Having now described various embodiments of the invention in detail as required by the patent statutes, those skilled in the art will recognize modifications and substitutions to the specific embodiments disclosed herein. Such modifications are within the scope and intent of the present invention as defined in the following claims.

## Claims

1. A method for selective sampling to assess texture using magnetic resonance (MR) comprising:
transmitting a first RF pulse with a first gradient chosen for first slice selection;
transmitting a second RF pulse with application of a second gradient chosen for slice selective refocusing in a region defined by an intersection of the first slice and a second slice;
encoding a specific k-value with a selected gradient pulse;
transmitting a third RF pulse with a third gradient activated, said third gradient adapted for slice selective refocusing, defining a region defined by the intersection of the first and second slices and a third slice selection to define a volume of interest (VOI);
turning off all gradients;
recording multiple samples of an RF signal encoded with the specific k-value;
applying a non-zero magnitude gradient as a time-dependent phase-encode determining a trajectory through k-space while recording samples at a sequence of k-values across a neighborhood of k-values defined by height and pulse width of the non-zero magnitude gradient, the sequence of k-values being a subset of k-values required to make an image;
post processing samples at a combination of sequential k values, recorded within a time span while the non-zero magnitude gradient is applied.

2. The method as defined in claim 1 wherein the step of encoding comprises:
applying an encoding gradient pulse to induce phase wrap to create a spatial encode for the specific k-value and orientation, the specific k-value selected based on texture within the VOI.

3. The method as defined in claim 2 wherein the specific k value is in a range of 0.2 mm−1 to 100 mm−1 and the neighborhood of k-values induced by the non-zero magnitude gradient are a subset of k-values required to make an image.

4. The method as defined in claim 2 further comprising:
applying a second non-zero magnitude gradient as a time-dependent phase-encode determining a second trajectory through k-space while recording samples at a second set of sequential k-values in a second neighborhood of the specific k-value, with all recording in a single TR.

5. The method as defined in claim 4 further comprising:
transmitting a second refocusing RF pulse prior to applying the second non-zero magnitude gradient.

6. The method as defined in claim 4 wherein the second non-zero magnitude gradient comprises a gradient inducing samples in k-space which are a subset of the k-values required to make an image.

7. The method as defined in claim 4 wherein the specific k value is in a range of 0.2 mm−1 to 100 mm−1 and the second non-zero magnitude gradient comprises a gradient inducing samples in k-space outside the neighborhood of k-values induced by the non-zero magnitude gradient, the k-values in the second neighborhood induced by the second non-zero gradient remain a subset of k-values required to make an image.

8. The method as defined in claim 2 further comprising:
following the recording of multiple samples of an RF signal encoded with the specific k-value, issuing a k-value selection pulse on a selected vector combination gradient to determine a second k-value;
turning off the vector combination gradient; and,
recording multiple samples of the RF signal at the second k-value in a single TR.

9. The method as defined in claim 8 further comprising:
applying additional k-value selection pulses for a predetermined plurality of pulses on a selected vector combination gradient, each k value selection pulse determining a next k-value;
turning off the vector combination gradient after each pulse; and,
recording multiple samples of the RF signal at the next k-value determined by each k value selection pulse in the TR.

10. The method as defined in claim 1 further comprising:
transmitting a refocusing RF pulse prior to applying the non-zero magnitude gradient.

11. The method as defined in claim 1 wherein the non-zero magnitude gradient is determined based on a windowing function and is defined as




  G
  ≤
  
    
      
        
          3
          !
        
        ⁢
        
          (
          
            1
            -
            
              η
              min
            
          
          )
        
      
    
    
      π
      ⁢
      
        
      
      ⁢
      
        X
        ·
        N
        ·
        Δ
      
      ⁢
      
        
      
      ⁢
      t
    
  




wherein ηmin is a lower correlation bound, X is a constant determining the windowing function and N is the number of samples.

12. The method as defined in claim 1 further comprising correcting for change in k-values across the neighborhood in k-space selected by the non-zero magnitude gradient by inducing a set of k-value measurements with prescribed k-value encodes.

13. A method for selective sampling to assess texture using magnetic resonance (MR) comprising:
transmitting RF pulses and associated gradients to define a volume of interest (VOI);
applying an encoding gradient pulse to induce phase wrap to create a spatial encode for a specific k-value and orientation determined based on texture within the VOI as an initial value;
applying a non-zero magnitude gradient as a time-dependent phase-encode determining a trajectory through k-space while recording samples at a sequence of k-values across a neighborhood proximate the specific k-value defined by height and pulse width of the non-zero magnitude gradient, the sequence of k-values being a subset of k-values required to make an image; and,
post processing the samples at the sequence of k values, recorded within a time span while the non-zero magnitude gradient is applied.

14. The method as defined in claim 13 wherein recording samples of a set at a sequence of k-values comprises sampling in overlapping blocks of k-values.

15. The method as defined in claim 14 wherein samples at the sequence of k-values in the blocks are combined as data points for an improvement in SNR via the combination and further comprising:
characterizing a power distribution across a range of the combined data points to measure underlying texture within the VOI.

16. The method as defined in claim 13 wherein the specific k value is in a range of 0.2 mm−1 to 100 mm−1 and the sequence of k-values in the neighborhood are a subset of the range.

17. The method as defined in claim 13 wherein the subset of k-values is less than 30% of k-values required to make an image.

18. A method for selective sampling to assess texture using magnetic resonance (MR) comprising:
exciting and internal volume of interest (VOI);
applying an encoding gradient pulse to induce phase wrap to create a spatial encode for a specific k-value and orientation, the specific k-value selected based on expected texture within the VOI;
applying a time varying gradient determining a trajectory through k space, from the specific k-value, of a sequence of k-values in a neighborhood of the specific k-value; and,
recording a set of signal measurements at the sequence of k-values in the neighborhood to measure texture within the VOI.