Famis: A software package for functional feature extraction from biomedical multidimensional images

0895-61 I l/92 $5.00 + .oO Copyright 0 1992 Pergamon Press Ltd.

Compurerized Medical Imaging and Graphics. Vol. 16. No. 2. pp. 81-91. 1992 Printed in the U.S.A. All rights reserved.

FAMIS: A SOFTWARE PACKAGE FOR FUNCTIONAL FEATURE EXTRACTION FROM BIOMEDICAL MULTIDIMENSIONAL IMAGES FredCrique Frouin

*f ,

Jean-Pierre Bazin*, Mireille Di Paola$, Odile Jolivet*, and Robert Di Paola*

*U66 INSERM, Institut Gustave-Roussy, Villejuif, France, and *Service de MCdecine NuclCaire, Institut Gustave-Roussy, Villejuif, France (Received 2 I June 1990; Revised 29 October 1991) Abstract-An increased number of image sequences is acquired in all modalities of the biomedical imaging field in order to study displacement or metabolism of a tracer or a contrast agent. It requires effective processing methods to estimate the underlying physiological components. We have developed a software package based on factor analysis algorithms which can adapt to various imaging modalities and its extension to double-indexed image sequences. We describe general characteristics of the software and present the main points of the user-friendly

interface. The performances of the package are discussed and the possibilities of the methodology are illustrated using an example in magnetic resonance imaging. Key Words: Functional imaging, Factor analysis, Image sequence, Medical imaging software, Multidimensional data

of interest (ROIs) method computes the average time activity or time concentration curve inside a functional structure previously delimited by the user. Thus the information is reduced to a few characteristic curves. ROIs are available in most of the biomedical image processing packages. Parametric imaging, in contrast with ROIs, aims at condensing information into a few images, which describe the local fit of the curve associated with one pixel to an a priori model (2), the retained parameters being presumably related to a physiological mechanism. The limitations of these two methods have been well studied. The lack of reproducibility of ROIs is due to both the interobserver and intraobserver variability in drawing regions. The alternative of automatic delimitation of ROIs does not provide better results and strongly depends on the acquisition protocol (3). Moreover, subtraction procedures are required to overcome the spatial structural overlapping. The major drawbacks of parametric imaging are its sensitivity to noise and its lack of universality, due to the adaptation of the model to each application (4). In conclusion, both methods reduce the information but not exhaustively. However, Fourier analysis of equilibrium radionuclide angiocardiography, which computes the amplitude and the phase of the first harmonic, is nearly exhaustive as long as the assumption of the sinusoidal model is valid (5). The method of Factor Analysis of Dynamic Structures (FADS) (6-8) is more general since it estimates the whole set of physiological components from any dynamic studies without making any modelling

INTRODUCTION

Increasing interest is presently devoted to functional imaging, which covers all the modalities of biomedical field. Recently Adam (1) defined three classes of tasks that can be assessed by functional imaging: organ motion, excretory functions, and metabolism steps. Except for organ motion, such studies require the injection of a tracer or a contrast agent. This domain was initiated in nuclear medicine and has now been extended to other modalities, for instance magnetic resonance imaging (MRI) using paramagnetic agents like chelates of gadolinium and computerized tomography (CT) employing iodide contrast agents. Time-dependent image sequences constitute the basic data set for functional imaging since they record both the spatial and dynamic distribution of the tracer or contrast agent. Three methods are mainly used to analyse these dynamic frames: visual inspection, regions of interest drawing, and parametric imaging. The simplest one, which is used in everyday clinical routine, consists of displaying the whole set of images either on an X-ray film or a video monitor, possibly using a cinematic mode. However, it does not provide quantitative and reproducible indexes. Computer-aided methods, aiming at a reliable reduction of information in the image sequence, were developed to overcome these drawbacks. The regions

+Correspondence should be addressed to Frederique Frouin, U66 INSERM, Institut Gustave-Roussy, 39, rue Camille Desmoulins, F-94805 Villejuif Cedex, France. 81

82

Computerized Medical Imaging and Graphics

March-April/

1992, Volume

16. Number 2

assumption about kinetics. The algorithm is based on an orthogonal transform, the Principal Component Analysis (PCA), followed by an oblique rotation under positivity constraints. It automatically computes factors, estimates of the underlying kinetics, and factor images that correspond to the spatial distributions associated with each factor. FADS is operator independent, exhaustive, and provides a good separation of the spatial structures even if they are superimposed. FADS has been implemented for nuclear medicine imaging on SIMIS S4000 minicomputers and distributed by Informatek Inc. since 1982. To date other groups working in the same domain are using FADS only for planar scintigraphy (9-13). Since 1982, we have proposed three ways to extend the method to other modalities:

where S(p, u) is the signal recorded at the pixel p for the value n of the observation variable (time or energy). K represents the number of underlying structures, Ck(u) is the kth kinetics or spectrum, I&) the spatial distribution of the corresponding structure, and E(p, v) an error term of the model. Positivity constraints are applied both to image and kinetics or spectral terms. The goal of the processing is to estimate the K functions G(u), called factors, and their associated images Zk(p), called factor images. Analogous monodimensional models are encountered in chemistry to estimate, for instance, the different components of mixtures by means of spectrometric measurements. A similar resolution approach based on PCA and oblique rotation under positivity conditions was proposed (17).

1) Processing of energy-dependent image sequences in scintigraphy ( 14): the method is similar to FADS, except that the time dimension is replaced by that of energy. It has recently been applied to multiple isotope immunoscintigraphic studies ( 15). 2) Processing of dynamic studies from various imaging modalities: CT, MRI, positron emission tomography (PET), fluorescence microscopy ( 16). 3) Analysis of multiple indexed image sequences ( 16) like: a) time- and energy-indexed image sequences in scintigraphy to process dynamic multi-isotope studies; b) time- and depth-indexed image sequences in MRI or PET to simultaneously study dynamic mechanisms in multiple slices.

Main steps of the algorithm

We have built a new software package, called FAMIS for Factor Analysis of Medical Image Sequences, to process these different types of data. It was implemented on general mainframes equipped with an image processor for display functions. This paper provides a complete description of the software. FACI-OR ANALYSIS OF MEDICAL IMAGE SEQUENCE Underlying model

We assume that the number of mechanisms present in the image sequence is reduced. Each one is characterized by a specific function, kinetics or spectrum, whose shape is identical for all points of the associated structure. As the spatial structures may overlap, the signal measured at each pixel is a weighted sum of the underlying kinetics or spectra. The formula (1) describes this basic hypothesis (8): S(P, V) = 2 Zk(P)- G(v) + mh k=l

VI

(1)

The algorithm proceeds in four steps (7, 8): data preprocessing, PCA, oblique rotation under positivity constraints, and factor image computation (see Fig. 1). Three-dimensional pixels, called “trixels,” are defined as the evolution of the signal according to the observation variable. Data preprocessing consists first of clustering trixels according to a rectangular pattern in the image, then of selecting the most significant ones by a threshold procedure based on their average values in the sequence. This procedure improves the signal to noise ratio and decreases the computing time. Selected trixels are scaled, centered, and submitted to a PCA. Prior scaling procedure implies that the K factors belong to the vector subspace generated by the (K 1) first principal components (7, 8). Each trixel is projected in this subspace and the K factors are estimated by an iterative oblique rotation under positivity constraints. The positivity conditions on factor images imply that the factors are the apices of the polytope that contains all the trixel projections. Moreover, positivity conditions on factor coefficients imply that the K factors belong to the convex set of the trixels having positive or null coordinates. This set always includes the origin (i.e., the average trixel). Once the K factors are estimated, the factor images are computed by projecting all the trixels of the original sampling on the oblique basis of the K factors. Double-indexed

image sequence

The extension of FAMIS to a global analysis of double-indexed image sequences was recently described (16). It requires that the underlying model can be expressed by the formula (2): S(P, VI 3 212)

=

S(q,

u2)

=

i k=l

($44)

- Cdv2)

+

@a

~2)

(2)

FAME

GENERAL

Selection of patient and data set 0 Preprocessing step : rectangular

Scaling and centring of data 0 Principal Component Analysis hk

.

Factor image (or image sequence) computhion

V I

FAMIS result display

I

End of processing

SOFI-WARE

FEATURES

First, the FAMIS software package handles various imaging modalities from different manufacturers, a situation encountered in many medical imaging departments. This implies dealing with a multiplicity of data storage formats and a variety of frame sizes. In addition, FAMIS works on both image sequences and double-indexed image sequences. Second, the software is portable. It was implemented on two different systems in our laboratory and distributed for external collaborations. To minimize the adaptation time of the different packages, we defined a general architecture for high-level image processing software (Fig. 2). Finally the package is flexible: it can run almost automatically in clinical applications and is at the same time highly interactive to solve research problems. Therefore a userfriendly interface was created. FORTRAN 77 was chosen as the programming language because of its extended use in the scientific community. The package was built around four lowlevel libraries devoted to mathematical functions, file input/output management, image and graphics display, user interface, and dialogue (Fig. 2). Each library is system-dependent, but the higher-level code is completely transportable. The proportion of FAMIS code that has to be modified for each new system is less than 5%.

aggregation

+I

83

F. FROWN et al.

??

I

Fig. 1. Flowchartof FAMIS software (general overview). where p is the space variable, 0, and v2 the two observation variables; q is a composite variable representing the space variable p and the variable vI . For instance, concerning time- and energy-indexed sequences in scintigraphy, v2 is the energy variable and q is a mixed variable including space and time dimensions. For time- and depth-indexed sequences (in MRI or PET), v2 is the time variable and q is a combination of the three space variables. v2 represents the main observation variable. Similarly to trixels, four-dimensional pixels, or “quadrixels,” are defined as the evolution according to the main variable 02 of the signal at a given value of q. Data preprocessing consists of clustering quadrixels according to a rectangular pattern and in thresholding resulting clusters. The two steps of PCA and oblique rotation are not modified. The final computation of Gk(q) leads to K factor image sequences that depend on the variable vI. In case of dynamic multislice sequence (in MRI or PET), FAMIS estimates K physiological kinetics Ck(v2) and associated factor images, Gk(q), computed for each slice.

Mathematical library The package uses some modules of a general mathematical library which is partitioned into two subsets: 1) the two-dimensional image-oriented library for elementary functions of image processing, and 2) the arithmetical library which includes floating point operations on vectors and matrices (e.g., inversion, singular value decomposition, etc.). These functions are generally programmed using standard instructions. However, if the systems are equipped with coprocessors or array processors, the mathematical procedures can be rewritten using the specificity of these devices. Input/output requirements Acquisition data (images and image sequences) are translated into a unique format, specific to the image or image sequence definition of each computer system. As an example, Table 1 shows the different sources already converted to perform FAMIS clinical applications. The lower-level functions of the data input/ output library define standard procedures to create, open, read, write, close, and delete files. Moreover, intermediate procedures were created to read/write elementary data sets like curves, images, and more complex files like image sequences or FAMIS results from

84

Computerized Medical Imaging and Graphics

March-April/ 1992, Volume 16, Number 2

Table 1. Current frame formats used for FAMIS applications

disk to memory or from memory to disk. A text file which summarizes the main parameters of the analysis and numerical results is also built, in order to be further edited, printed, or archived.

Imaging equipment Nuclear medicine

General Electric Starcam 500 A Sopha Medical DSX LET1 (PET)

CT

General Electric 9800 Elscint Exe1 2400

MRI

General Electric Signa CGR Magniscan Philips Gyroxan

Display and graphics procedures

A library of two-dimensional graphics and image display was created to unify the graphical libraries proposed by the manufacturers. We defined a set of elementary primitives which verify the following principles:

1) Each function is associated with a virtual window

curves or the display of a set of frames. These functions are completely device-independent since they were generated from primitives of the previously defined library.

of 4096 X 4096 pixels. 2) Only absolute coordinates are taken into account. Relative coordinates and the concept of current position do not exist. 3) The graphical elements are lines, connected segments of lines, portions of ellipses, filled areas, symbols, texts, and images. Each element is characterized by a set of attributes as is its position in the window, its size, color, or orientation.

User interface

The user interface library manages three kinds of interaction: 1) Values are entered with the keyboard. The fit of numerical values within lower and upper limits is verified while the alphanumerical values are filtered, to match the expected ones. 2) Scrolling menus were created to select the appropriate processing stages. The implementation used the specific libraries of screen management of each system. 3) The keyboard, which is the default device, interacts with the video monitor and, as an option, a specialized interface such as a mouse or a track ball.

The monitor is considered as a virtual window which can be divided into square or rectangular windows. Each square window can be identically split. Procedures were defined to assign the monitor to the current process, to select a window prior to any display order, and to choose associated parameters: look-up tables (LUT) for image and graphics. Higherlevel procedures were created to achieve all the graphical functions of FAMIS, for instance the drawing of

HIGH-LEVEL APPLICATION

r

0 0

I

100%compatible Specific of the material cotiguration

-

Necessary Optioru32

Fig. 2. General architecture of FAMIS package.

FAME

??

85

F. FROUIN et al.

Fig. 3. MRI osteosarcoma double-indexed image sequence after Gd-DTPA injection. Twelve dynamic frames were acquired in 4 slices: Sl, S2, S3, S4. The first frame is obtained before the injection. Memory management

aim of the exam is to assess the response

For new applications like dynamic MRI, CT, or double-indexed image sequences, both computing time and memory allocation are critical. A global optimization is performed according to the size of memory of each system and the characteristics of the current study. To do so, the arrays are dimensioned according to the main parameters of the processing, that is, size of frames, number of frames and modality number of each index for double-indexed image sequences, maximum number of factors, and number of regions after preprocessing. According to each configuration, indexes are computed to estimate whether the two largest arrays of the processing step (initial image sequence and factor images) can be stored in main memory. Otherwise, the image sequence must be read twice from the disk, and factor images should be temporarily stored on disk.

coma under chemotherapy. The 256 X 256 frames were acquired during 10 min after gadolinium-diethylenetriaminopenta-acetic acid (Gd-DTPA) injection, using a TZ-weighted sequence that simultaneously records the signal into four nearly contiguous 5 mm slices. The final set of frames was a double-indexed image sequence with 12 time intervals and 4 depth intervals (Fig. 3). For each slice, the first image, acquired just before the injection of Gd-DTPA is subtracted from the following ones in order to keep only the information related to the kinetics of the contrast agent.

USER INTERACTIONS

An example of dynamic multislice MR images illustrates the interactions provided by the package. The

of osteosar-

Data preprocessing

First the operator chooses the parameters (height and length) of the rectangular pattern of the spatial clustering. Then the cluster selection is achieved using the “sum image,” defined as the sum of all the images of the sequence. For double-indexed image sequence, V, sum images are defined, one for each value of the variable ul. Figure 4 shows the four sum images and the rectangular pattern of 8 X 8 pixels of the MRI

Computerized Medical Imaging and Graphics

March_April/l992,

Volume 16, Number 2

of external constraints. The first possibility is the display of any active trixel. It can be chosen either on the sum image or in the reduced subspace by pointing the coordinates of its projection on each axis. Its positions in both sum image and reduced subspace are then highlighted and the corresponding kinetics or spectrum is displayed in the third quadrant. Moreover, the function associated with any point of the reduced subspace can be displayed, too. The last option allows definition of constraints for the oblique rotation by delimiting the center and the radius of a “hypersphere” in the reduced subspace. This possibility was used to introduce external constraints as a feedback of compartmental analysis ( 18). Factor and factor image display Fig. 4. FAME processing of the frames of Fig. 3. The first frame of each slice is subtracted from the following ones. Final display of the four sum images (one per slice) and selected clusters after rectangular clustering of 8 X 8 pixels and thresholding procedure.

Figure 7 shows the three factors and their factor image sequences, one image per slice, that are suitable to describe our example. The three functional structures are interpreted as follows (19):

1) a peak followed by a fast decreasing component study. A cross indicates the suppressed clusters. The automatic mode is a thresholding procedure. The manual selection helps its refinement: the operator can switch the current status of any cluster by positioning the cursor on it and validating a control key. Oblique rotation

Two execution modes are available to perform the iterative step of oblique rotation: the automatic mode, in which the criteria that stop the iterations are a priori defined by using either default values or user defined values, and the manual mode, in which the operator decides whether to go on or to stop at each iteration. To facilitate the user control of the procedure, the current estimation of the factors is displayed. The subspace generated by the first (K - 1) principal components is represented by a set of planes, each one being related to a couple of principal components (Fig. 5). The projections of the trixels are drawn in these different planes. The edges of the polytope, the apices of which are the current factor estimates, are displayed, too. Exploration of the reduced subspace

This option, which is highly interactive, connects the spatial position of any trixel in the sum image(s) with its position in the reduced subspace and displays the corresponding kinetics or spectrum. The monitor is divided into four quadrants (Fig. 6): 1) the current sum image, 2) the reduced subspace, 3) the display of the functions associated to trixels, and 4) the display

which indicates an early arrival of the contrast agent in arteries and viable tumor; 2) a slower kinetics with a later peak followed by a decreasing component, which corresponds to late vascularization; 3) a kinetics of accumulation showing the diffusion of the tracer in extracellular interstitial tissues, which reveals the extension of oedematous and inflammatory mechanisms. Factor image contribution corresponds to the percentage of total information associated with the structure. Finally the superimposition of up to three different factor images may be displayed, using the three fundamental colors (i.e., red, green, and blue). The intensity of each component may be modulated according to the value of its contribution. Figure 8 synthesizes the whole information of the multislice MRI study employing a three-color superimposition of the three factor images, repeated for the four slices. This display emphasizes the overlapping of the spatial structures in spite of the thickness of the slices and demonstrates the ability of FAMIS to achieve the separation of superimposed structures. FAMIS IMPLEMENTATION Material configuration

The FAMIS package was implemented on three systems composed of the following elements: 1) a VAX/VMS

8300 (Digital Equipment Corporation); an array processor MP32 (Floating Point Systems); and a SIGMEX 7000 graphics processor,

FAMIS

F. FROUIN et al.

87

??

which can manage up to 1024 X 1024 2%bit pixels, and a 1024 X 1024 monitor, using either pseudocolor with 256 colors from a palette of 16 million or true color with 8 bits per fundamental. In addition, four planes are dedicated to graphical overlay. 2) a MV4OOO/AOS (Data General); an array processor MP32 (Floating Point Systems); an image processor (General Electric), which can manage up to eight planes of 5 12 X 5 12 X 16 bits, connected to a video processor of 5 12 X 5 12 pixels with 8 bits for the image display, allowing pseudocolor or monochromatic displays and 4 bits for graphical overlay. For two-color or three-color display, a LUT was implemented where each fundamental is encoded using from 6 to 16 levels. 3) a VAX750/VMS (Digital Equipment Corporation); and a graphics processor (ADAGE), the characteristics of which are very similar to the GE image processor. Performance evaluation The performances of the algorithm were evaluated on VAX 8300. Central processing unit (CPU) times are given in Table 2. The tested configurations correspond to conventional biomedical applications. As our computer system was equipped with an array processor, the computation of factor images was implemented using the mathematical library provided by the manufacturer. Table 3 enables us to compare elapsed times measured in identical conditions to compute factor images with and without using the array processor. Elapsed times, which correspond to the operator’s real times spent at the console, are considerably decreased for large sequences. DISCUSSION

What we have learnedfrom our experience with FAMIS Numerous and various clinical applications of FAMIS were undertaken using either the FADS soft-

E

external curves

curve *

trixe11 w

-1

Fig. 6. Schematics of the monitor for exploration of reduced subspace; upper left: sum image and selected trixels, upper right: position of trixels, factors, and constraints in the subspace generated by the first principal components, lower left: display of trixels previously selected, lower right: display of external constraints introduced in the analysis.

ware on SIMIS computers for scintigraphic data (2027) or the new FAMIS package for dynamic (multislice) MRI (19), dynamic CT (28), dynamic fluorescent microscopy ( 16), and scintigraphic image sequences depending on time and energy (29). These studies have proved the superiority of results obtained by FADS to those estimated by ROI’s methods. Our approach has taught us that a learning stage is necessary for each new clinical protocol in order to define the optimal parameters of analysis, particularly the number of factors, the adequate size of the rectangular clustering, and the threshold applicable to trixels. Therefore, several analyses must be performed on the PC24

limits of the positive

PC31

--6$/p*’ PC2

I

K=3 factors

: proJeCtIon Ot trixels PCi : ith principal component

Fi

*

+

: estimation of the ith factor.

Fig. 5. Display of oblique rotation.

PC3

I K=4 factors

Computerized Medical Imaging and Graphics

March-April/ 1992, Volume 16, Number 2

same study training set to test different values of the parameters. Moreover, this practice is useful in learning how to interpret FAMIS results and to study how normal and pathological situations are revealed on factors and factor images. For instance, the interpretation of FAME results in the case of MRI osteosarcoma studies was assessed by the comparison of factor images with histological maps of resected specimen ( 19). Once the processing protocol is set up, the sofiware can be used very simply in a quasi-automatic mode. Furthermore, the interactive exploration of the reduced subspace is a powerful tool in understanding the limits of the method. It allows us to test solutions to overcome these limits, since a priori known kinetics or spectra may be introduced in this reduced subspace. They can be manually imposed like real factors (16, 18) to improve estimation. At the present time we are working on refining algorithms by introducing a priori information. This new approach is currently recommended (9-l 3) to overcome some difficulties inherent in the conventional method. Improvement of the user interface The package needs to be developed further to improve file management, graphical output, and its integration into a heterogeneous network. In our environment the main limitation of FAMIS use is the overloading of graphics devices. We intend to direct display tasks toward small workstations offering a

Fig. 7. Factors and factor image sequences estimated from the MRI sequence by FAME. Four factor images corresponding to the four slices are associated with each factor, Fi. The three factors correspond to early vascular (Fl), late vascular (F2) and interstitial kinetics (F3).

Fig. 8. Three-color superimposition of spatial structures estimated by FAME. Red: early vascular, blue: late vascular. green: interstitial.

FAMIS

F. FROWN et al.

89

??

Table 2. VAX 8300 CPU times for FAMIS execution Numberof framesV,. V, Modality Nuclear medicine Nuclear medicine MRI, CT Nuclear medicine MRI

Size of frames

V,

v2

Number of clusters

Number of factors

CPU time (min:s)

64 X 64 128 x 128 256 x 256 64 X 64 256 X 256

1 1 1 20 4

60 60 11 20 11

100 400 1,000 2,000 4,000

4 4 3 4 3

1:45 6:30 4:40 22:05 22:15

V2 = number of modalities of the principal observation variable; VI = number of modalities of the second observation variable (time or slice) (V, is equal to 1 when there is only one observation variable).

ing. For instance, dynamic frames computed from a scintigraphic sequence depending both on time and energy can be further analyzed by FAMIS (29). Data storage as well as the development of image databases must take into account this particularity of processed results to successfully integrate source images and their processing. Part of these improvements have already been integrated in the new version of FAMIS for clinical systems, currently being developed by Sopha Medical, a nuclear medicine manufacturer.

minimum of display functions and to reserve the highquality graphics processor for three-color superimposition. The development of inexpensive graphics boards should solve this problem. The second limitation is associated with data communication and file management. Up to now most of data are first collected on acquisition devices, then loaded onto processing machines via magnetic tapes. The development of Picture Archiving and Communication Systems (PACS) inside and outside hospitals should improve data exchange. Concerning file management, we are working on a physical data storage independent of higher level applications. At the present time, a logical file which specifies the path to physical data has been implemented only for MRI data. It is derived from the characteristics of the acquisition sequence given by the operator. In the future, the information of the file identifiers will be managed by a database to create these logical files automatically (30). The management of FAMIS results is more complex. The whole set of results: K factors (curves) and K factor images or factor image sequences form an entity, which is related to the original sequence and can also be used for display and/ or processing. Moreover each element or group of elements (e.g., curve, frame, or image sequence) may be submitted to a further processing like filtering, curve fitting, quantification, image comparison, or process-

SUMMARY Functional imaging, which is now available in various biomedical imaging modalities (e.g., nuclear medicine, MRI, and CT), can be considered as an indispensable complement of anatomical and morphological observations to assess a diagnosis. It currently uses time-dependent image sequences to follow up the excretion and/or the metabolism of a tracer or a contrast agent. Few methods analyze these image sequences and extract reproducible and quantitative indexes. Conventional techniques are based on either time concentration curve extraction inside a region of interest or parametric image computation. FAMIS cumulates the goals of both methods because it yields the curves

Table 3. Elapsed times measured to compute and store factor images with and without using FPS MP32 array processor

Nuclear medicine Nuclear medicine MRI, CT Nuclear medicine MRI

Size of frames

V,

64 X 64 128 X 128 256 X 256 64 X 64 256 X 256

1 1 20 4

1

v2

60 60 11 20 11

Number of factors

Frames in memory

Elaspsed time on VAX 8300

Elapsed time on FPS

4 4 3 4 3

Yes No No No No

0:32 3:30 5:30 12:05 17:30

0:26 0:45 0:30 3:50 1:45

V2 = number of modalities of the principal observation variable; V, = number of modalities of the second observation variable (time or slice) I when there is only one observation variable).

(V, is equal to

90

Computerized Medical Imaging and Graphics

and the associated images that estimate the underlying physiological mechanisms. In contrast to ROIs, the method is largely automatic and exhaustive. Compared to parametric imaging, it is less noise sensitive and it does not require any assumption about the shape of underlying elements. The algorithm is based on a PCA followed by an oblique rotation under positivity constraints. The resulting software package fits the initial description of the method: FADS to process dynamic radionuclide studies. The package was developed to process multiple imaging modalities. Furthermore, the FAMIS methodology was improved to analyze doubleindexed image sequences like time and energy sequences in nuclear medicine or multislice dynamic sequences in MRI. The architecture of the software package can be extended to any high-level application of biomedical imaging. It is based on device-dependent libraries concerning mathematical, display, user interface, and file management functions. Moreover, it provides a synthetic display at each step of processing, which allows an interactive modification of main parameters of analysis. The software is designed to run quasi-automatically in clinical applications and to allow high interactivity in research, in order to investigate and overcome some intrinsic limitations of the method. Most of the clinical applications were initially performed in scintigraphy. However, FAMIS interest in other imaging modalities has been recently demonstrated, like the example of dynamic multislice sequences in MRI to study the response of osteosarcoma to chemotherapy. The computing times of the algorithm were evaluated on VAX 8300. They are already realistic in the case of dynamic sequences. They are reasonable concerning double-indexed image sequences and could be easily decreased using current technologies. However, such studies require handling a large amount of data. In the future, data exchange should be facilitated using efficient networks and databases to manage both original and processed images.

Acknowledgments-This work was supported by the Association pour la Recherche sur le Cancer (ARC) and Sopha Medical Inc. We wish to thank Dr. A. Chamentier for supnlvina MRI studv. We are erateful . to 1. Kuchenthal for-her translation help.

REFERENCES 1. Adam, W.E. A general comparison of functional imaging in nuclear medicine with other modalities. Semin. Nucl. Med. 17:317; 1987. 2. Goris, M.L. Parametric images as a tool for quantitative normative evaluation. Semin. Nucl. Med. 17: 18-27; 1987.

March-April/ 1992, Volume 16, Number 2 3. Cosgriff, P.S. Region of interest, ROI, analysis-Confidence in derived parameters? Nucl. Med. Comm. 6:305-309: 1985. 4. Holman, B.L. Merits and limitations of functional imaging techniques. In: Schmidt, H.A.E.; Rosier, H., eds. Nuklearmedizin: Computer assisted functional analysis. Stuttgart, Germany: F. K. Schattauer Verlag; 1982: 9-l 3. 5. Bacharach, S.L.; Green, M.V.; Vitale, D.; White, G.; Douglas, M.; Bonow, R.; Larson, S. Optimum Fourier filtering of cardiac data: A minimum-error method. J. Nucl. Med. 26:1176-l 184; 1983. 6. Bazin, J.P.; Di Paola, R.; Gibaud, B.; Rougier, P.; Tubiana, M. Factor analysis of dynamic scintigraphic data as a modelling method: An application to the detection of metastases. In: Di Paola, R.; Kahn, E., eds. Information processing in medical imaging, vol. 88. Paris: lnstitut National de la Sonti et de la Recherche Medicale (INSERM); 1980:345-366. I. Barber, D.C. The use of principal components in the quantitative analysis of gamma camera dynamic studies. Phys. Med. Biol. 25~283-292; 1980. 8. Di Paola, R.; Bazin, J.P.; Aubry, F.; Aurengo, A.; Cavailloles, F.; Herry, J.Y.; Kahn, E. Handling of dynamic sequences in nuclear medicine. IEEE Trans. Nucl. Sci. NS29: 13lO- 1321: 1982. 9. Nijran, K.S.; Barber, DC. Factor analysis of dynamic function studies, using a priori physiological information. Phys. Med. Biol. 31:1107-1117; 1986. 10. Nijran, K.S.; Barber, DC. The importance of constraints in factor analysis of dynamic studies. In: De Graaf, C.; Vietgever, M., eds. Information processing in medical imaging. New York: Plenum Press; 1988:521-529. 11. Houston, AS. The use of cluster analysis and constrained op timisation techniques in factor analysis of dynamic structures. In: Mathematics and computer science in medical imaging. Berlin: Springer Verlag; N.A.T.O. A.S.I. Series F, 1988:39:491-503. 12. Samal, M.; Kamy, M.; Surova, H.; Marikova, E.; Dienstbier, Z. Rotation to simple structure in factor analysis of dynamic radionuclide studies. Phys. Med. Biol. 32:371-382; 1987. 13. Nakamura, M.; Suzuki, Y.; Kobayashi, S. A method for recovering physiological components from dynamic radionuclide images using the maximum entropy principle: A numerical investigation. IEEE Trans. Biomed. Eng. 36(9):906-917; 1989. 14. Di Paola, R.; Bazin, J.P.; Aubry, F.; Kahn, E.; de Talhouet, H.; Aurengo, A., Cavailloles, F.: Herrv, J.Y.: Kahn. E. Developments in data processing in nuclear medicine imaging. Izotoptechnika 25:94-l 1 I; 1982. 15. Manil, L.; Boudet, F.; Motte, P.; Gardet, P.; Saccavini, J.C.; Lumbroso, J.; Schlumberger, M.; Caillou, B.; Bazin, J.P.; Ricard, M.; Bellet, D.; Bohuon, C.; Tubiana, M.; Di Paola, R.; Parmentier, C. Positive anti-calcitonin immunoscintigmphy in patients with medullarv thyroid carcinoma. Cancer Research 49:54805485; 1989. - 16. Di Paola, R.; Frouin, F.; Bazin, J.P. Utilisation des sequences d’images dans les etudes metaboliques. In: Venot, A.; Degoulet, P., eds. Informatique et medicaments. Paris: Springer-Verlag; 1989:97-109. 17. Lawton, W.H.; Sylvestre, E.A. Self-modeling curve resolution. Technometrics 13:617-633; 197 1. 18. Frouin, F.; Bazin, J.P.; Di Paola, R. Image sequence processing using factor analysis and compartmental modelling. In: Vietgever, M., ed. Science and engineering of medical imaging. SPIE, 1137: 3.7-44; 1989. 19. Charpentier, A.; Bonnerot, V.; Frouin, F.; Vane], D.; Charlotte, F.; Kalifa, C.; Di Paola, R. Factor analysis processing of dynamic MR imaging: New method to assess osteosarcoma preoperative chemotherapy response. Radiology 177:22 1; 1990 (Abst.). 20. Cavailloles, F.; Bazin, J.P.; Di Paola, R. Factor Analysis in gated cardiac studies. J. Nucl. Med. 25:1067-1079; 1984. 2 1. Cavailloles, F.; Bazin, J.P.; Capderou, A.; Valette, H.; Hebert, J.L.; Di Paola, R. Automated processing of first-pass radionuclide angiocardiography by factor analysis of dynamic structures. Nucl. Med. Comm. 8:375-387; 1987. 22. Pavel, D.G.; Sychra, J.; Olea, E.; Kahn, C.; Virupannavars, S.; Zolnierczyk, K.; Shanes, J. Factor analysis and its place in the

FAMIS

23.

24. 25.

26.

27.

28.

29.

30.

??

F. FROWN et al.

evaluation of ventricular regional wall motion abnormalities. In: Bacharach, S.L., ed. Information processing in medical imaging. Dordrecht, The Netherlands: Martinus Nijho@ 1986: 193-206. Villanueva-Meyer, J.; Philippe, L.; Cordero, S.; Mardes, C.S.; Mena, I. Use of factor analysis in the evaluation of left to right cardiac shunts. J. Nucl. Med. 27:1442-1448; 1986. Philippe, L.; Mena, I.; Darcourt, J.; French, W.J. Evaluation of valvular regurgitation by factor analysis of first pass angiography. J. Nucl. Med. 29:159-167; 1988. Aurengo, A.; Bazin, J.P.; Lumbroso, J.; Di Paola R. Factor analysis as an aid in the interpretation of dynamic scintigrams of transplanted kidneys. In: Raynaud, C., ed. Nuclear medicine and biology advances, vol. 4. Elmsford, NY: Pergamon Press; 1982:3 167-3 170. Herry, J.Y.; Bazin, J.P.; Bourguet, P.; Di Paola, M.; Di Paola, R. Factor analysis of dynamic scintigraphic structures: Application to kinetics of hepatobiliary tracers. In: Raynaud, C., ed. Nuclear medicine and biology advances, vol. 2. Elmsford, NY; Pergamon Press; 1982:2240-2243. Cinotti, L.; Bazin, J.P.; Meignan, M.; Aubry, F.; Di Paola, R. Factor analysis of *lmKr lung ventilation studies. Eur. J. Nucl. Med. IO:51I-518; 1985. Chabriais J.; Di Paola, R.; Bazin, J.P.; Rougier, P. Towards blood flow quantization in hepatic masses by dynamic CT. Proceedings of the 17th International Congress of Radiology, Paris, 467; 1989. Helal, B.O.; Bazin, J.P.; Goeau, 0.; Leguillouzic, D.; Archambaud, F.; Lebtaki, R.; Pueyo, M.; Desgrez, A.; Di Paola, R. Diagnostic de la malignite des nodules is&s non fonctionnels au moyen d’une Ctude dynamique en double marquage. J. de MCdecine NuclCaire et Biophysique 13:79-96; 1989. Aubry, F.; Kaplan, H.; Di Paola, R. An image handling system for medical image processing. In: Viergever, M., ed. Science and engineering of medical imaging. SPIE, 1137:131-138; 1989.

About the Author-FReDfiRIQUE FROUINgraduated from the Ecole

Nationale Sup&ieure des T&communications (ENST), Paris, in 1986. She received her Ph.D. in Signal Processing from ENST in 1989. From 1986 to 1989 she attended the Institut de Formation Sup&ieure BiomCdicale, Villejuif, France, to follow a theoretical

91

training in medical and biological engineering. In 1990, she was a post-doctoral fellow at Sopha Medical Inc., But, France, where she implemented FAMIS soAware on nuclear medicine workstations. She is currently Char& of Research of the Institut National de la Santt et de la Recherche MCdicale (INSERM) in the Unit& 66. Her main interests are biomedical image processing and physiological modelling. About the Author-JEAN-PIERRE BAZINworked at the Commissariat B I’Energie Atomique (CEA) from 1961 to 1968, where he was involved in biological system modelling. He has been Research Engineer of INSERM in the Unit& 66 since 1968. He has also worked both on compartmental modelling and on image sequence processing. He has been particularly concerned with the development of factor analysis techniques applied to image sequence processing. About the Author-MIREILLE DI PAOLAhas been working as a computer scientist in the Nuclear Medicine Department of the Institut Gustave-Roussy since 1975. She has been involved in several project realizations concerning the development of nuclear medicine workstations. She is currently working on multimodality medical image processing and picture archiving solutions. About the Author-ODILE

JOLIVETis Agrkgte in chemistry and she received the D.E.A. degree in computer science. She has been Engineer at the Centre National de la Recherche Scientifique (CNRS) since 1978. She joined the Unit& 66 of INSERM in 1985. She is involved in computer science and image processing, particularly in the field of MRI.

About the Author-ROBERT DI PAOLA was born in Toulouse, France, on May 4, 1940. He received the D.E.S. in physics from the University of Toulouse in 196 1. He joined the Nuclear Medicine Department at the Institut Gustave-Roussy, Villejuif, France, in 1963, as Attach6 of Research of the Faculty of Medicine of Paris. In 1966, he became Attach6 of Research of INSERM in the Unit& 66. In 1969 he became ChargC of Research and in 1985 he became Director of Research in the Unit& 66. Since 1986 he has led this Unit& which deals with morphological and functional biomedical imaging. He has specialized in the processing of medical images, both two-dimensional and threedimensional, with emphasis on functional imaging.