Image processing method, program, and recording medium

Abstract

The image processing method of the present invention includes: a step of performing optical coherence tomography on a cell block and acquiring three-dimensional image data representing a three-dimensional image of the cell block (step S111); a step of acquiring a thickness distribution representing the thickness of the cell block between the outer surface facing the external space and the inner surface facing the cavity in each position of the cell block based on the three-dimensional image data (step S116); and a step of detecting a protrusion whose thickness is greater than that of the surrounding area based on the acquisition result of the thickness distribution (step S118). Based on the OCT image data of a cell block having a cavity inside, quantitative information useful for analyzing the surface structure of the cell block can be obtained.

Description

Technical Field

[0001] The present invention relates to image processing for photographing a cell mass, particularly a cell mass having a cavity therein, and analyzing its structure. For example, the present invention relates to an image processing technique for providing information useful to a user evaluating a fertilized egg based on data corresponding to a three-dimensional image obtained by tomographically photographing the fertilized egg. Background Art

[0002] For example, in reproductive assisted medicine for the purpose of infertility treatment, embryos (fertilized eggs) that have been fertilized in vitro and cultured for a certain period of time are returned to the body. However, the pregnancy success rate in reproductive assisted medicine is not necessarily high, and the mental and economic burden on patients is also large. In order to solve this problem, a method for accurately determining the status of cultured embryos is explored.

[0003] For example, in an embryo that has developed to the blastocyst stage, a hollow structure called the "blastocyst cavity" is generated inside, and a layer of cells called the "trophectoderm" covers the periphery of the hollow structure. The number of cells constituting the trophectoderm is an indicator for evaluating the state of the embryo, and therefore, a technology for non-invasively counting the number of cells is required.

[0004] In the past, about the evaluation of whether embryo culture progresses well, usually by for example doctor, embryo culture personnel visually carried out by microscopic observation.As its judging index, for example Veeck classification, Gardner classification etc. are widely used.But these are only the rough judging benchmarks of expression for the morphological characteristics of embryo, and current situation is that final evaluation depends on the subjective judgment of assessor.Therefore, need can be objective and the technology of quantitative evaluation.

[0005] In addition, for example, it is known that the intestinal epithelial organoids are composed of a group of cells with cavities inside. When used as an evaluation model for a cell membrane transport system, it is ideal to form a cavity covered by a layer of cell layers. However, in reality, multiple layers of cell layers are often formed around the cavity. Therefore, it is necessary to quantitatively evaluate the structure of the cell layer on the surface of the intestinal epithelial organoids, specifically, how cells form a cell layer covering the cavity.

[0006] As a technology expected to meet this requirement, patent documents 1 and 2 previously applied for by the applicant of the present application are disclosed. In Japanese Patent Publication No. 2019-133429 (Patent Document 1), a method for identifying and segmenting the nutritional ectoderm and the inner cell mass from a three-dimensional image of an embryo (fertilized egg) photographed using non-invasive tomographic photography techniques such as optical coherence tomography (OCT) is recorded. In addition, in Japanese Patent Publication No. 2019-132710 (Patent Document 2), a technology for visualizing the embryo photographed by OCT as a two-dimensional image is recorded. Specifically, the following technology is recorded: the embryo is represented by polar coordinates with the center of gravity of the embryo as the origin, and the reflected light intensity in each radial direction is replaced by a brightness value to assist the observer in performing the evaluation operation of the embryo.

[0007] The above-mentioned prior art is a technology that presents the three-dimensional structure of the observed object in a manner that is easy for the observer to understand. However, it has not yet reached the point where effective quantitative information can be automatically extracted to assist in observation and evaluation operations. In this regard, there is room for improvement in the above-mentioned prior art. More specifically, there is a need for a technology that can obtain quantitative information related to the surface structure of a cell mass with a cavity inside, such as an embryo in the blastocyst stage and an intestinal epithelial organoid.

[0008] In the embryo (fertilized egg) of the blastocyst stage, as the main component of the evaluation object, there is also a zona pellucida, in addition to the above-mentioned nutrient ectoderm and inner cell mass. However, in patent documents 1 and 2, there is no mention of the technology of extracting the area corresponding to the zona pellucida. As a practical problem, in the embryo, the nutrient ectoderm is distributed in a manner pasted on the inner side of the zona pellucida, and there is almost no difference in the brightness information of the two in the OCT image. Therefore, it is difficult to distinguish the nutrient ectoderm and the zona pellucida from the OCT image.

[0009] As described above, there is no established method for clearly distinguishing the vegetative ectoderm from the zona pellucida in OCT images, and the techniques described in the above-mentioned Patent Documents 1 and 2 do not address this point. Therefore, in order to provide useful information to a user who wants to evaluate them separately, a technique is required that can more accurately distinguish the vegetative ectoderm and the zona pellucida in OCT images. Summary of the invention

[0010] The present invention is proposed in view of the above problems, and its purpose is to provide a technology that can obtain quantitative information useful for analyzing the surface structure of a cell block based on image data obtained by OCT in image processing of an image of a cell block with a cavity inside. For example, its purpose is to identify the vegetative ectoderm and the zona pellucida by using the image data of the embryo obtained by OCT, thereby effectively assisting the user in performing the evaluation operation of the embryo.

[0011] In order to achieve the above-mentioned purpose, an image processing method according to one embodiment of the present invention takes an image of a cell block having a cavity inside as a processing object, and includes: a process of performing optical coherence tomography on the cell block and acquiring three-dimensional image data representing a three-dimensional image of the cell block; a process of acquiring a thickness distribution representing the thickness of the cell block between an outer surface facing an external space and an inner surface facing the cavity in each position of the cell block based on the three-dimensional image data; and a process of detecting a protrusion having a thickness protruding compared to a surrounding area based on a result of acquiring the thickness distribution.

[0012] Due to the internal structure of the cell such as the nucleus, each cell has a structure in which the central part is expanded compared to the peripheral part. Therefore, in the cell layer formed by the cells being connected to each other laterally, bumps and depressions are generated due to the shape of each cell. Especially in a cell layer composed of a layer of cells in a manner covering the internal cavity, as in an embryo (fertilized egg) in the blastocyst stage, protrusions corresponding to each cell appear in the layer. Therefore, the number of protrusions can be an indicator of the position and number of cells constituting the layer. In addition, the intervals between the protrusions can be an indicator of the size of each cell.

[0013] In the above invention, the distribution of the thickness of the cell block between the outer surface of the cell block and the inner surface of the cavity facing the inside is obtained based on the three-dimensional image data of the cell block obtained by OCT. This thickness distribution can represent the thickness of the cell layer covering the cavity at each position. The thickness of the cell layer does not necessarily have to be uniform. For example, in a cell layer composed of a layer of cells, the central part of each cell is thick and the peripheral part is thinner. Therefore, if a protrusion with a larger thickness than the surrounding area is detected based on the obtained thickness distribution, it can be said that the possibility that this part is a part representing the position of a cell is high.

[0014] As described above, by obtaining the thickness distribution of the cell layer covering the internal cavity and detecting the protrusions included therein, it is possible to obtain information quantitatively indicating the position, number, size, etc. of the cells constituting the cell layer.

[0015] In addition, in order to achieve the above-mentioned purpose, an image processing method of one embodiment of the present invention includes: a process of obtaining three-dimensional image data obtained by optical coherence tomography of a fertilized egg in the blastocyst stage and two-dimensional image data obtained by optical microscopy; a process of determining a structural region corresponding to the structure of the fertilized egg based on the three-dimensional image data; a process of determining a region corresponding to the zona pellucida of the fertilized egg based on the two-dimensional image data and calculating its average thickness; and a process of dividing the structural region into a region within the range of the average thickness from the outer edge toward the inside of the structural region and other regions.

[0016] Detailed description as described later, especially, the outer shape of the embryo (fertilized egg) in the blastocyst stage of mammal is roughly spherical, and its surface is covered by the zona pellucida with roughly constant thickness.And, at the inner surface of the zona pellucida, the nutritional ectoderm with various thicknesses is distributed in the mode of pasting according to the position.In addition, in the three-dimensional image of the embryo obtained by optical coherence tomography (OCT), the difference of brightness is hardly seen between the nutritional ectoderm and the zona pellucida.That is, in the OCT image, it is only difficult to distinguish the nutritional ectoderm and the zona pellucida by brightness value.

[0017] Therefore, in the present invention, by combining the use of OCT images and optical microscope images, it is possible to distinguish between the nutritional ectoderm and the zona pellucida in the OCT image. Specifically, based on the OCT image, the structural region corresponding to the structure of the embryo (fertilized egg) is determined. At this time, there is no need to distinguish between the nutritional ectoderm, the zona pellucida, and the inner cell mass, as long as the region occupied by a certain structure in the three-dimensional image of the fertilized egg is determined. Here, if based on the insight that the thickness of the zona pellucida is roughly constant, it can be assumed that the range of a certain distance from its outer edge toward the inside in the structural region determined by the OCT image is occupied by the zona pellucida.

[0018] On the other hand, the thickness of the zona pellucida is estimated from an image obtained by photographing the same fertilized egg with an optical microscope. In an optical microscope image, it is easier to identify the vegetative ectoderm and the zona pellucida based on the difference in brightness between the vegetative ectoderm and the zona pellucida than in an OCT image. Taking advantage of this, the thickness of the zona pellucida can be estimated from the optical microscope image.

[0019] Based on this information, the nutritional ectoderm and the zona pellucida can be distinguished in the OCT image. That is, in the three-dimensional image of the fertilized egg obtained by OCT shooting, the range corresponding to the thickness of the zona pellucida obtained from the optical microscope image from its outer edge can be regarded as the area occupied by the zona pellucida. On the other hand, in the structure that cannot be distinguished from the zona pellucida in the OCT image, the area that is more inner than the range corresponding to the thickness of the zona pellucida from the outer edge can be regarded as the area occupied by the nutritional ectoderm.

[0020] As described above, by dividing the three-dimensional image of the fertilized egg obtained by OCT into a region within the range of the thickness of the zona pellucida derived from the optical microscope image from the outer edge toward the inside of the structure region and other regions, the zona pellucida and other structures can be distinguished.

[0021] Another embodiment of the present invention is a program for causing a computer to execute the above-mentioned steps. Another embodiment of the present invention is a computer-readable recording medium that non-transitorily stores the above-mentioned program.

[0022] As described above, according to the present invention, quantitative information related to the position and number of each cell, etc., useful for analyzing the surface structure of the cell mass can be obtained based on the image data obtained by performing OCT imaging of a cell mass having a cavity inside. For example, by combining the information obtained from the OCT image with the information obtained from the optical microscope image, the area occupied by the zona pellucida and other areas can be identified with high accuracy in the OCT image of the embryo (fertilized egg). BRIEF DESCRIPTION OF THE DRAWINGS

[0023] Figure 1 This is a diagram showing a configuration example of an image processing device suitable as an execution subject of the image processing method of the present invention.

[0024] Figure 2A This is a diagram schematically showing the structure of an embryo serving as a sample in this embodiment.

[0025] Figure 2B This is a diagram schematically showing the structure of an embryo serving as a sample in this embodiment.

[0026] Figure 3 : is a flowchart showing the image processing in this embodiment.

[0027] Figure 4 This is a flowchart showing an example of a specific method for constructing a classification model.

[0028] Figure 5A It is a graph showing the relationship between the focal depth and the area of ​​the extracted transparent zone.

[0029] Figure 5B It is a graph showing the relationship between the focal depth and the area of ​​the extracted transparent zone.

[0030] Figure 6 : is a flowchart showing a process for selecting an image focused on a transparent zone.

[0031] Fig. 7A It is a diagram schematically showing how the region is divided according to the present embodiment.

[0032] Figure 7BIt is a diagram schematically showing how the region is divided according to the present embodiment.

[0033] Figure 8 This is a diagram showing the correspondence between the XYZ orthogonal coordinate system and the polar coordinate system.

[0034] Fig.9A This is a diagram illustrating a method for calculating the thickness of the trophectoderm.

[0035] Fig. 9B This is a diagram illustrating a method for calculating the thickness of the trophectoderm.

[0036] Fig.10 This is a flowchart showing the process of obtaining the thickness of the trophectoderm.

[0037] Fig.11A This is a diagram showing an example of a two-dimensional mapping method.

[0038] Fig. 11B This is a diagram showing an example of a two-dimensional mapping method.

[0039] Fig. 12A This is a diagram showing another example of the two-dimensional mapping method.

[0040] Fig. 12B This is a diagram showing another example of the two-dimensional mapping method.

[0041] Fig.13A is a diagram showing an example of a three-dimensional mapping method.

[0042] Fig. 13B is a diagram showing an example of a three-dimensional mapping method.

[0043] Fig.14A It is a diagram showing the principle of peak detection processing.

[0044] Fig. 14B It is a diagram showing the principle of peak detection processing.

[0045] Fig. 14C It is a diagram showing the principle of peak detection processing.

[0046] Fig.15 is a flowchart showing the peak detection process.

[0047] Fig.16A : is a diagram showing an example of the peak position detected on the two-dimensional graph.

[0048] Fig. 16B : is a diagram showing an example of the peak position detected on the two-dimensional graph.

[0049] Fig.17A : is a diagram showing an example of the peak position detected on the three-dimensional graph.

[0050] Fig. 17B : is a diagram showing an example of the peak position detected on the three-dimensional graph. DETAILED DESCRIPTION

[0051] Figure 1 1 is a diagram showing the structure of an image processing device suitable as an execution subject of the image processing method of the present invention. The image processing device 1 performs tomographic photography of a sample carried in a liquid, for example, an embryo (fertilized egg) cultured in a culture medium. Then, the obtained tomographic image is image processed, and a cross-sectional image showing the structure of a cross section of the sample is produced. In addition, a stereoscopic image of the sample is produced based on a plurality of cross-sectional images. In order to indicate the direction of the drawings, as shown in FIG. Figure 1 As shown in FIG. 1 , an XYZ orthogonal coordinate axis is set. Here, the XY plane represents a horizontal plane. In addition, the Z axis represents a vertical axis, and more specifically, the (-Z) direction represents a vertically downward direction.

[0052] The image processing device 1 includes a holding unit 10. The holding unit 10 holds a sample container 11 containing a sample S as a photographic object in a horizontal position. The sample container 11 is, for example, a container in the shape of a flat-bottomed shallow dish called a "dish". Alternatively, for example, the sample container 11 is a well plate having a plurality of recesses (holes) capable of carrying a liquid formed on the upper surface of a plate-like member. A culture medium M such as a culture solution is injected into the sample container 11, and the sample S is carried inside the culture medium M.

[0053] An imaging unit 20 is provided below the sample container 11 held by the holding unit 10. An optical coherence tomography (OCT) device capable of capturing a tomographic image of an object in a non-contact, non-destructive (non-invasive) manner is used in the imaging unit 20. As described in detail later, the OCT device, i.e., the imaging unit 20, includes a light source 21 for generating illumination light for an object, a fiber coupler 22, an object optical system 23, a reference optical system 24, a spectrometer 25, and a photodetector 26.

[0054] The photographing section 20 also includes a microscope photographing unit 28 for performing optical microscope photography. More specifically, the microscope photographing unit 28 includes a photographing optical system 281 and a photographing element 282. The photographing optical system 281 includes an objective lens, and the focus of the objective lens is focused on the sample S in the sample container 11. In addition, as the photographing element 282, for example, a CCD photographing element, a CMOS sensor, etc. can be used. As the microscope photographing unit 28, a microscope photographing unit capable of bright field photography or phase difference photography is preferably used. The object optical system 23 and the microscope photographing unit 28 are supported by a support member (not shown) that can be moved in the horizontal direction, and can change the position in the horizontal direction.

[0055] The image processing device 1 further includes a control unit 30 for controlling the operation of the device and a driving unit (not shown) for driving a movable portion of the imaging unit 20. The control unit 30 includes a CPU (Central Processing Unit) 31, an A / D converter 32, a signal processing unit 33, an imaging control unit 34, an interface (IF) unit 35, an image memory 36, and a memory 37.

[0056] The CPU 31 controls the operation of the entire device by executing a predetermined control program. The control program executed by the CPU 31 and the data generated during the processing are stored in the memory 37. The A / D converter 32 converts the signal outputted from the photodetector 26 of the imaging unit 20 and the imaging element 282 according to the amount of received light into digital data. The signal processing unit 33 performs the signal processing described below based on the digital data outputted from the A / D converter 32, and produces an image of the object. The image data produced in this way is appropriately stored and stored by the image memory 36.

[0057] The photographing control unit 34 controls each part of the photographing unit 20 so that it performs photographing processing. Specifically, the photographing control unit 34 selectively positions the object optical system 23 for OCT photographing and the microscope photographing unit 28 for optical microscope photographing at a specified photographing position that contains the photographed object, i.e., the sample S, in the photographing field of view. When the object optical system 23 is positioned at the photographing position, the photographing control unit 34 causes the photographing unit 20 to perform the OCT photographing processing described later, and obtains three-dimensional image data representing the three-dimensional structure of the sample S. On the other hand, when the microscope photographing unit 28 is positioned at the photographing position, the photographing control unit 34 causes the microscope photographing unit 28 to obtain two-dimensional image data representing a planar image of the sample S formed on the light receiving surface of the photographing element 282 through the photographing optical system 281.

[0058] The interface unit 35 is responsible for the communication between the image processing device 1 and the outside. Specifically, the interface unit 35 has a communication function for communicating with external devices and a user interface function for accepting input operations from the user or notifying the user of various information. In order to achieve this purpose, an input device 351 and a display unit 352 are connected to the interface unit 35. The input device 351 is composed of a keyboard, a mouse, a touch panel, etc. that can accept input operations related to the function selection, action condition setting, etc. of the device. In addition, the display unit 352 is composed of a liquid crystal display, etc. that displays various processing results such as tomographic images and stereoscopic images produced by the signal processing unit 33.

[0059] In the imaging unit 20, a low-coherence light beam including a broadband wavelength component is emitted from a light source 21 having a light-emitting element such as a light-emitting diode or a superluminescent diode (SLD). For the purpose of imaging a sample such as a cell, for example, near-infrared light is preferably used so that the incident light reaches the inside of the sample.

[0060] The light source 21 is connected to an optical fiber 221, which is one of the optical fibers constituting the optical fiber coupler 22. The low-coherence light emitted from the light source 21 is branched by the optical fiber coupler 22 into light to two optical fibers 222 and 224. The optical fiber 222 constitutes the object system optical path. More specifically, the light emitted from the end of the optical fiber 222 is incident on the object optical system 23.

[0061] The object optical system 23 includes a collimating lens 231 and an objective lens 232. The light emitted from the end of the optical fiber 222 is incident on the objective lens 232 via the collimating lens 231. The objective lens 232 has the function of converging the light (observation light) from the light source 21 on the sample and the function of converging the reflected light emitted from the sample and directing it toward the fiber coupler 22. A single objective lens 232 is shown in the accompanying drawings, but it can also be a combination of multiple optical elements. The reflected light from the object is incident on the optical fiber 222 as signal light via the objective lens 232 and the collimating lens 231. The optical axis of the objective lens 232 is orthogonal to the bottom surface of the sample container 11. In this example, the direction of the optical axis is consistent with the direction of the vertical axis.

[0062] CPU31 issues a control command to the shooting control unit 34. In response, the shooting control unit 34 moves the shooting unit 20 in a specified direction. More specifically, the shooting control unit 34 moves the shooting unit 20 in the horizontal direction (XY direction) and the vertical direction (Z direction). By moving the shooting unit 20 in the horizontal direction, the shooting range changes in the horizontal direction. In addition, by moving the shooting unit 20 in the vertical direction, the focal position of the objective lens 232 in the optical axis direction changes relative to the object to be photographed, i.e., the sample S.

[0063] A portion of the light incident on the fiber coupler 22 from the light source 21 is incident on the reference optical system 24 via the optical fiber 224. The reference optical system 24 includes a collimating lens 241 and a reference mirror 243. These constitute a reference system optical path together with the optical fiber 224. Specifically, the light emitted from the end of the optical fiber 224 is incident on the reference mirror 243 via the collimating lens 241. The light reflected by the reference mirror 243 is incident on the optical fiber 224 as reference light.

[0064] The reference mirror 243 is supported by an advance and retreat member (not shown) that is actuated by a control command from the imaging control unit 34, and can be moved forward and backward in the Y direction. The reference mirror 243 moves in the Y direction, that is, in a direction approaching and moving away from the collimating lens 241, thereby adjusting the optical path length of the reference light reflected by the reference mirror 243.

[0065] The reflected light (signal light) reflected by the surface or internal reflective surface of the sample and the reference light reflected by the reference mirror 243 are mixed in the fiber coupler 22 and incident on the photodetector 26 via the optical fiber 226. At this time, interference caused by the phase difference occurs between the signal light and the reference light, but the spectral spectrum of the interference light is different depending on the depth of the reflection surface. That is, the spectral spectrum of the interference light has information in the depth direction of the object. Therefore, by splitting the interference light for each wavelength and detecting the light amount, and Fourier transforming the detected interference signal, the reflected light intensity distribution in the depth direction of the object can be obtained. The OCT shooting technology based on this principle is called "Fourier Domain OCT (FD-OCT)".

[0066] In the imaging unit 20 of this embodiment, a spectrometer 25 is provided on the optical path of the interference light from the optical fiber 226 to the photodetector 26. As the spectrometer 25, for example, a spectrometer using a prism, a spectrometer using a diffraction grating, etc. can be used. The interference light is split into wavelength components by the spectrometer 25 and received by the photodetector 26.

[0067] Based on the interference light detected by the light detector 26, the interference signal output from the light detector 26 is Fourier transformed, thereby obtaining the reflected light intensity distribution in the depth direction of the incident position of the illumination light in the sample, that is, in the Z direction. By scanning the light beam incident on the sample container 11 in the X direction, the reflected light intensity distribution on the plane parallel to the XZ plane is obtained. Based on the result, a tomographic image of the sample with the plane as a cross section can be produced. Since this principle is well known, a detailed description is omitted.

[0068] In addition, by changing the incident position of the light beam in the Y direction in multiple stages while taking tomographic images at each position, a large number of tomographic images of the sample taken in cross sections parallel to the XZ plane can be obtained. If the scanning pitch in the Y direction is reduced, image data with a resolution sufficient to grasp the three-dimensional structure of the sample can be obtained. Based on these tomographic image data, three-dimensional image data (so-called voxel data) corresponding to the three-dimensional image of the sample can be produced.

[0069] Thus, the image processing device 1 has a function of acquiring an image of the sample S carried in the sample container 11 together with the culture medium M. As images, two-dimensional image data acquired by an optical microscope, tomographic image data acquired by OCT, and three-dimensional image data based on tomographic image data can be acquired.

[0070] Hereinafter, an embodiment of the image processing method of the present invention that can be executed using the image processing device 1 configured as described above is described. The image processing of this embodiment has the content of extracting the region corresponding to the main structures constituting the embryo, namely, the zona pellucida, the vegetative ectoderm, and the inner cell mass, from an image of a fertilized egg (hereinafter, simply referred to as an "embryo") at the blastocyst stage taken as a sample S.

[0071] Based on the data thus obtained, it is possible to effectively assist users (specifically physicians or embryo culture personnel) in performing embryo evaluation operations. For example, in the culture of embryos for the purpose of infertility treatment, the image processing method of this embodiment can be applied for the purpose of obtaining insights for determining whether the culture is well performed.

[0072] Figure 2A and Figure 2B It is a diagram schematically showing the structure of an embryo used as a sample in this embodiment. As is known to all, when an egg is fertilized, cleavage begins and the blastocyst is formed through a state called "morula". Figure 2A Schematically shows the internal structure of an embryo in the blastocyst stage. In the blastocyst stage, the embryo E has a cavity called "blastocyst cavity B" inside. More specifically, cells undergoing cleavage are arranged as a thin layer on the surface of the embryo and form a vegetative ectoderm T, and the internal space surrounded by the vegetative ectoderm T forms the blastocyst cavity B. In addition, an inner cell mass I in which a large number of cells are densely packed is formed in a part of the internal space.

[0073] In addition, the zona pellucida ZP is formed so as to cover the outer surface of the vegetative ectoderm T. The zona pellucida ZP is a membrane with a substantially uniform thickness mainly composed of glycoprotein. On the other hand, the vegetative ectoderm T formed by a large number of cell aggregations has various thicknesses depending on the position and is distributed so as to adhere to the inner surface of the entire zona pellucida ZP. Figure 2B As shown in FIG. 1 , the vegetative ectoderm T is formed by a large number of cells C arranged along the inner surface of the zona pellucida ZP. Since the shapes and sizes of the individual cells C are different, the vegetative ectoderm T becomes a layer with different thicknesses depending on the position. Figure 2A In the present study, the deviation in the thickness of the trophic ectoderm T was emphasized more than the actual value.

[0074] These structures that constitute embryo E are zona pellucida ZP, trophectoderm T and inner cell mass I as the area of ​​interest in the evaluation of embryo and are of great significance.Therefore, the technology of automatically extracting the zone corresponding with these structures from the image photographed is of great significance in assisting the user to carry out the evaluation operation of embryo.Yet, in the tomographic image or the three-dimensional image obtained by OCT shooting, the brightness information between these structures almost does not have difference.Therefore, in the separation carried out according to simple brightness difference, it is difficult to separate these structures with high accuracy.

[0075] Thus, the technology for automatically separating the zona pellucida ZP, the vegetative ectoderm T, and the inner cell mass I from the OCT image has not been established so far. In view of such a problem, the applicant of the present application previously disclosed patent document 2 that can separate the regions corresponding to the vegetative ectoderm and the inner cell mass from the OCT image. However, even if this technology is adopted, it cannot be said that the zona pellucida is clearly separated from the vegetative ectoderm. Specifically, the region extracted as the vegetative ectoderm by this technology may include the region corresponding to the zona pellucida.

[0076] In the image processing of this embodiment, the OCT image is analyzed using information obtained from the two-dimensional image obtained by taking an optical microscope, thereby separating the zona pellucida ZP and the vegetative ectoderm T. Then, the first purpose is to separately determine the areas occupied by the zona pellucida ZP, the vegetative ectoderm T, and the inner cell mass I in the three-dimensional image. In addition, the second purpose is to calculate the number of cells constituting the vegetative ectoderm T separated in this way.

[0077] In the bright field image or phase contrast image obtained by optical microscope shooting, it is relatively easy to distinguish zona pellucida ZP and vegetative ectoderm T according to the difference of its brightness and texture. Then, in the embryo of the blastocyst stage well cultivated, zona pellucida ZP has a roughly uniform thickness. In other words, it can be considered that in the embryo, the region from its surface to a certain depth is occupied by zona pellucida ZP. Thus, the thickness of zona pellucida ZP can be determined according to the two-dimensional image data obtained by optical microscope shooting. Then, according to the three-dimensional image data obtained by OCT shooting, it is determined that the region corresponding to zona pellucida ZP can be separated and extracted from vegetative ectoderm T from the surface observation of embryo E. In addition, for the separation of inner cell mass I, for example, the method of recording in patent documentation 2 can be applied. Through these, the first purpose can be achieved.

[0078] like Figure 2B As shown, the vegetative ectoderm T separated from other structures in this way has unevenness caused by the shape of each cell C. More specifically, protrusions caused by the shape of each cell C are generated on the surface of the vegetative ectoderm T facing the blastocyst cavity B. Therefore, if such protrusions are detected and the number is counted, the number of cells C constituting the vegetative ectoderm T can be known. In this way, the second purpose can be achieved.

[0079] Figure 31 is a flowchart showing the image processing of the present embodiment. The processing is realized by the CPU 31 executing a control program prepared in advance so that each part of the device performs a predetermined action. When the sample container 11 containing the embryo as the evaluation object is placed in the holding part 10 (step S101), the embryo is taken as the object to be photographed, and optical microscope photography and OCT photography performed by the photography part 20 are performed.

[0080] exist Figure 3 In the embodiment, steps S102 to S105 represent the optical microscope shooting and the data processing using the two-dimensional image data obtained thereby. On the other hand, steps S111 and S112 represent the OCT shooting and the image processing using the three-dimensional image data obtained thereby. These two processes can be performed one after the other in time series, in which case it is arbitrary which process is performed first. In addition, it is also possible to shorten the processing time by performing these processes in parallel.

[0081] First, the summary of the processing of steps S102 to S105 is described. In this processing, the two-dimensional image data of the embryo is obtained by photographing the embryo with an optical microscope, and the thickness of the zona pellucida is calculated based on the image data. Specifically, the microscope photographing unit 28 is positioned to the shooting position, and the focus position is changed in multiple stages in the depth direction (Z direction), and photographed each time. Thus, a collection of multiple two-dimensional images with different focus depths, i.e., a so-called Z stack image (step S102), is obtained.

[0082] The region corresponding to the zona pellucida ZP is extracted from each of these images (step S103). Then, among the multiple images, the image that is best focused on the zona pellucida ZP is selected (step S104). The thickness of the zona pellucida ZP is calculated based on the selected image (step S105). The following will describe these processes in detail, that is, the process of extracting the region corresponding to the zona pellucida ZP (step S103), the process of selecting the most focused image (step S104), and the process of calculating the thickness of the zona pellucida ZP (step S105). Through the processes so far, the thickness of the zona pellucida ZP in the embryo E of the evaluation object is known.

[0083] On the other hand, in step S111, the object optical system 23 is positioned at the shooting position, and the embryo E is tomographically photographed while scanning the shooting position. Thus, the three-dimensional image data of the embryo E is obtained. The brightness of each part represented by the three-dimensional image data is binarized with a specified threshold value (step S112). Thus, the three-dimensional image of the embryo E is separated into a high-brightness area (hereinafter referred to as "structure area") occupied by some structures of higher density and a low-brightness area of ​​lower density. For example, the culture medium M is low-brightness in the OCT image.

[0084] As described above, the structures occupying the structural region may include the zona pellucida ZP, the vegetative ectoderm T, and the inner cell mass I. In steps S113 and S114, these structures are separated from each other. Specifically, in step S113, the region corresponding to the zona pellucida ZP is separated from the other structural regions. Then, in step S114, the vegetative ectoderm T is separated from the inner cell mass I.

[0085] Through the processing so far (steps S101 to S114), the first purpose of determining the areas occupied by the zona pellucida ZP, the vegetative ectoderm T, and the inner cell mass I from the three-dimensional image of the embryo E is achieved. The processing of the following steps S115 to S118 is for the second purpose of counting the number of cells constituting the separated vegetative ectoderm T. The processing for this purpose will be described in detail later.

[0086] The following describes the various element technologies for executing the various processes (steps S103 to S105, S111) of the above-mentioned processing in sequence. It should be noted that since the process of obtaining the Z stack image of step S102 is well known, its description is omitted. In addition, the process of step S114 can also be applied to the technology described in patent document 2, so its description is omitted.

[0087] In step S103, extract the area corresponding to the zona pellucida ZP from the two-dimensional image data that the embryo E is taken with an optical microscope. Suitable image processing techniques can be used to perform this processing. For example, the pattern recognition technology in the area with determined features can be applied. Specifically, the classification model is constructed by using the image of the zona pellucida obtained in advance as the learning with a teacher of the teacher image, and the optical microscope image of the embryo E of the object of evaluation is carried out to regional segmentation with this. By this, the area corresponding to the zona pellucida ZP can be extracted from the image.

[0088] As a specific example of the region segmentation process, for example, a well-known semantic segmentation method can be used. The semantic segmentation method is a technique for labeling each pixel in an image using a classification model pre-built by a deep learning algorithm. In the present embodiment, this method can be used as follows.

[0089] First, the user prepares the optical microscope image of the embryo that the zona pellucida is photographed with good image quality, and each pixel in the region corresponding to the zona pellucida in the image is labeled to represent the meaning. Then, by using the original microscope image as input data and the label image as the correct solution data to perform deep learning, a classification model is constructed. If an unknown image is provided as input data to the classification model pre-built in this way, an output image of the label of the meaning in the region corresponding to the zona pellucida in the image can be obtained. By extracting such a region from the output image, the zona pellucida can be extracted as a result.

[0090] Figure 4 1 is a flowchart showing an example of a specific method for constructing a classification model. The processing can be performed by various computer devices having a function of displaying images and a function of accepting input operations from users. For example, it can be performed by a general-purpose computer device such as an image processing device 1 or a personal computer.

[0091] First, display the optical microscope image (step S201) of the embryo that is pre-photographed in the state of focusing (focusing) on ​​the transparent zone. In the image processing device 1, it is possible to make the display unit 352 display the image. For the image displayed in this way, accept the teaching input (step S202) from the user for specifying the area suitable for the transparent zone of focusing. In this case, it is expected that the user is a skilled person who has enough insights into the image of the embryo. In addition, when using the image processing device 1, it is possible to accept the teaching input via the input device 351.

[0092] The region designated as the transparent zone is given a label indicating its meaning (step S203). The image labeled in this way is used as the correct answer data, and the original image is used as the input data to perform deep learning. In this way, a classification model for extracting the transparent zone from the image is constructed (step S204). Labels other than the transparent zone can also be used as needed.

[0093] This classification model is a model constructed by using an image focused on the transparent band as an input image. Therefore, when it is applied to an unknown test image and the semantic segmentation method is performed, an area with strong features as a transparent band is extracted from the test image. When the transparent band is extracted in this way, its thickness can be calculated. In order to calculate the thickness with high accuracy, it is desired that as large an area as possible within the image is focused on the transparent band. In other words, it can be said that the calculation of the thickness of the transparent band is desired to be performed based on an image in which the area extracted by the semantic segmentation method is as large as possible. It should be noted that the area can be represented by the number of pixels, for example.

[0094] On the other hand, the embryo E of the evaluation object has a three-dimensional structure. Therefore, in the image that is appropriately determined to have a focal depth and is taken, there may be a situation where the focusing state of the zona pellucida ZP is not necessarily good. Therefore, from the Z stack images obtained at different focal positions in the depth direction (Z direction), an image with the largest area of ​​the region extracted as a region equivalent to the zona pellucida ZP is selected. Then, the image is used for the calculation of the thickness of the zona pellucida ZP.

[0095] In step S104, an image that is most focused on the zona pellucida ZP is selected from the Z stack image. In the semantic segmentation method of this embodiment, the region of the zona pellucida in focus is extracted from the image. Therefore, it can be said that the image with the largest area in this region is likely to be the image that is most focused on the zona pellucida ZP. An image that meets such conditions can be selected.

[0096] Figure 5A as well as Figure 5B is a graph showing the relationship between the focal depth and the area of ​​the extracted transparent zone. Figure 5A The upper part of schematically shows the relationship between the optical microscope images (Z stack images) taken at each focal position in the depth direction and the regions (regions corresponding to the zona pellucida ZP) extracted from these images by the semantic segmentation method. In addition, Figure 5A The curve graph at the bottom of illustrates the relationship between the focal position during shooting and the area (number of pixels) of the extracted region. As shown in these, in the best-focused image Ia in the Z stack image with a clear transparent zone ZP, the area of ​​the extracted region is the largest. In other words, in the Z stack image, the focal position corresponding to the image with the largest area of ​​the extracted region can be regarded as the focal position. Therefore, this image Ia can be selected as the basic image for calculating the thickness of the transparent zone ZP.

[0097] However, blurring may occur in the image due to vibration during shooting, for example, and as a result, the transparent zone may appear larger in apparent area than in reality. Figure 5B 2 is a diagram showing such a situation. For example, in image Ib, image blurring due to vibration occurs. As a result, the area extracted as a region corresponding to the transparent zone becomes larger in appearance, and this position may be misjudged as the focus position.

[0098] In this embodiment, this problem can be solved by using the brightness difference between pixels of the peripheral portion of the extracted area. That is, it is considered that in a well-focused image, the boundary between the area corresponding to the transparent band and the surrounding area is clear, so there is also a clear contrast between the brightness of these areas. On the other hand, in an unfocused image, the boundary between these areas is not clear, so the contrast is not sharp.

[0099] Thus, instead of evaluating only based on the size of the area (number of pixels) of the extracted region, an evaluation value that also reflects the size of the brightness change of the edge portion of the region, i.e., the sharpness, is introduced. By doing so, it is believed that the above-mentioned misjudgment can be reduced. There are many methods for quantifying such edge changes, and these methods can be appropriately selected and applied.

[0100] In this embodiment, as an example, the value obtained by multiplying the area of ​​the extracted region by a coefficient reflecting the magnitude of the brightness change of the edge portion is used as the evaluation value. As the coefficient, for example, a value obtained by squaring the brightness difference between pixels across the edge can be used. More specifically, the difference between the average brightness of all pixels adjacent to the edge in the extracted region and the average brightness of all pixels adjacent to the edge outside the edge is calculated, and the coefficient is calculated by squaring the difference.

[0101] By doing so, Figure 5B As shown in the lower curve, the peak is more distinct, and the risk of the increase in the extraction area due to the blur in the image Ib becoming a cause of misjudgment of the focus position can be reduced. It should be noted that the brightness difference is squared here to make the coefficient a positive value, but instead, the absolute value of the brightness difference can also be used as the coefficient.

[0102] Figure 6 is a flowchart showing a process for selecting an image focused on a transparent zone, the content of which is related to Figure 3 The processing equivalent to step S104 is performed. For each image constituting the Z stack image, the area equivalent to the transparent zone is extracted in step S103. In order to find the area of ​​the area extracted in this way, the number of pixels belonging to the area is counted for each image (step S301). Then, the average brightness of the adjacent pixels inside and outside the edge of the extracted area is found, and the difference is calculated (step S302). Based on these values, the evaluation value representing the focus of each image is calculated (step S303). Specifically, the evaluation value is calculated by multiplying the number of pixels in the extracted area by a coefficient represented by the square of the brightness difference inside and outside the edge. The image with the largest evaluation value found in this way is selected as the image that is most focused on the transparent zone.

[0103] Next, the processing content of step S105 for calculating the thickness of the zona pellucida ZP according to an optical microscope image selected from the Z stack image is described. For the selected image, the extraction of the area equivalent to the zona pellucida ZP is performed in step S103. In the well-cultured embryo E, it is considered that the annular area with a roughly constant width is extracted as the area equivalent to the zona pellucida ZP. In the image focused on the zona pellucida ZP, the width of the annulus, that is, the distance between the inner edge of the annulus and the outer edge corresponds to the thickness of the zona pellucida ZP.

[0104] Various methods can be considered as methods for finding the width of the ring. As a simple method, for example, there is a method of using the distance transform (Distance Transform) function equipped in the open source computer vision (OpenCV, Open Source Computer Vision) library. By taking a pixel on the inner edge of the ring as a focus pixel and applying the Distance Transform function, the distance from the pixel to the pixel on the outer edge of the nearest ring can be determined. This distance represents the width of the ring at this position, that is, the thickness of the transparent band ZP. Conversely, it is also equivalent to taking the pixel on the outer edge of the ring as a focus pixel and finding the shortest distance from the pixel to the inner edge. The average value of the width found at each position on the ring in this way can be used as the average value of the thickness of the transparent band ZP. In the following, the average value of the thickness is represented by the symbol Tav.

[0105] The above has described the detailed processing contents of steps S103 to S105. Next, the separation processing of the zona pellucida in step S113 will be described. Here, in the three-dimensional image of the embryo E obtained by OCT shooting, the range from the surface to the depth Tav is regarded as the zona pellucida ZP. Therefore, by extracting only this range from the three-dimensional image, it is possible to take out only the structure corresponding to the zona pellucida ZP. On the other hand, by removing the structure of this range from the three-dimensional image, it is possible to take out the structure outside the zona pellucida ZP, that is, take out the trophectoderm T and the inner cell mass I.

[0106] For the three-dimensional image excluding the region corresponding to the zona pellucida ZP, a region segmentation process using, for example, a Local Thickness operation is performed on the remaining structure, as described in Patent Document 2. In this way, the region corresponding to the vegetative ectoderm T and the region corresponding to the inner cell mass I can be separated. In this way, the regions occupied by the zona pellucida ZP, the vegetative ectoderm T, and the inner cell mass I can be extracted from the three-dimensional image of the embryo E, respectively.

[0107] Fig. 7A and Figure 7B Schematically shows the regional segmentation of this embodiment. Here, the embryo E is shown as a two-dimensional image for ease of understanding. However, in fact, it is a process for the three-dimensional image obtained by OCT. Fig. 7A As shown, the region ranging from the surface Es of the image of the embryo E to the thickness Tav is regarded as a region corresponding to the zona pellucida ZP.

[0108] In addition, if Figure 7BAs shown, when the local thickness (Local Thickness) operation is performed to the region suitable for removing the zona pellucida ZP from the three-dimensional image of embryo E, the region represented as the aggregate of smaller spheres is regarded as the vegetative ectoderm T, and the region represented as the aggregate of larger spheres is regarded as the inner cell mass I. This distinction can be realized as described below. That is, a suitable threshold value is set to the radius of the inscribed sphere generated by the Local Thickness operation. Then, the region represented by the ball with a radius less than this threshold value is regarded as the vegetative ectoderm T. On the other hand, the region represented by the ball with a radius larger than the threshold value is regarded as the inner cell mass I.

[0109] When it is necessary to further subdivide each region separated in this manner, as described in Patent Document 2, region segmentation can be performed on each separated region by, for example, a watershed method.

[0110] In this way, by using the information obtained by the three-dimensional image data obtained by OCT shooting of the embryo (fertilized egg) in the blastocyst stage and the two-dimensional image data obtained by optical microscope shooting of the same embryo, the zona pellucida in the structure constituting the embryo can be separated from other structures, especially from the nutritional ectoderm that is difficult to produce a difference in brightness with the zona pellucida in the OCT image. Specifically, the average thickness of the zona pellucida is calculated based on the optical microscope image that is prone to the difference in texture between the zona pellucida and the nutritional ectoderm, and the area within the range of the average thickness from its surface in the structure in the three-dimensional image is regarded as the area occupied by the zona pellucida. Using the situation that the thickness of the zona pellucida is roughly the same, the information of the thickness of the zona pellucida obtained from the optical microscope image is reflected in the three-dimensional image to perform the segmentation of the region. In this way, the zona pellucida and the nutritional ectoderm can be separated with high precision.

[0111] Next, the process for achieving the second object of the present embodiment, that is, counting the number of cells constituting the trophectoderm T cells, is described. Figure 3 The general flow of this process is to calculate the thickness distribution of each position of the trophectoderm T, detect the peak position in the obtained thickness distribution, and count the number of obvious peaks as peaks corresponding to each cell. The following describes each process in turn.

[0112] In step S115, the three-dimensional image data represented by the XYZ orthogonal coordinate system is converted into polar coordinates. Since the embryo E is roughly spherical and hollow inside, the vegetative ectoderm T has a shape close to a spherical shell. In order to more simply represent such a structure, it is preferred to represent each position with polar coordinates (spherical coordinates) with the center of the embryo E as the origin. Therefore, as Figure 8As shown, coordinate conversion is performed from the XYZ orthogonal coordinate system to the rθφ polar coordinate system having one radius vector r and two deflection angles θ and φ as coordinate variables.

[0113] Figure 8 It is a diagram showing the correspondence between the XYZ orthogonal coordinate system and the polar coordinate system. As we all know, if the coordinates (x, y, z) of the origin O at point P in the orthogonal coordinate system and the coordinates (r, θ, φ) of point P in the spherical coordinate system are common, then the following relationship exists:

[0114] x=r·sinθ·cosφ

[0115] y=r·sinθ·sinφ

[0116] z=r·cosθ

[0117] Specifically, the center position of embryo E is determined according to the three-dimensional image data, and the center position is used as the origin O of the spherical coordinate system. The origin O need not be consistent with the origin of the XYZ orthogonal coordinate system in the raw signal data. Then, the coordinate conversion from the orthogonal coordinate to the spherical coordinate is performed by appropriate conversion processing. By performing coordinate conversion in this way, each position in the three-dimensional space determined by the XYZ orthogonal coordinate system can be expressed by polar coordinates.

[0118] About " center of embryo ", based on three-dimensional image data, for example, can be obtained as described below. When the surface of the three-dimensional image of embryo E represented by three-dimensional image data is regarded as a spherical surface, the center of gravity of the ball in the image can be used as the center of embryo. In three-dimensional image processing, the method of obtaining the center of gravity of a solid object is well known, and this method can be applied. In addition, for example, it is also possible to determine a spherical surface or a rotating ellipsoidal surface that approximately represents the surface of the embryo, and use the center of the approximate surface as the center of the embryo.

[0119] In step S116, based on the three-dimensional image data converted into polar coordinates in this way, the thickness distribution of the trophectoderm T is obtained. Specifically, by calculating the thickness of the trophectoderm T in one radial direction in the polar coordinate space and performing the calculation for each radial direction, the thickness distribution of the trophectoderm T can be obtained.

[0120] Fig.9A and Fig. 9B This is a diagram showing a method for calculating the thickness of the trophectoderm. Fig.9A In the example shown, consider a radius vector r1 extending from the origin O. The distance T1 between the inner surface of the vegetative ectoderm T in the direction of the radius vector r1, that is, the surface Si facing the blastocyst cavity B and the outer surface, that is, the surface So facing the zona pellucida ZP, is taken as the thickness of the vegetative ectoderm T in the direction of the radius vector.

[0121] In addition, Fig. 9B In the example shown, the point P1 where the radius vector r1 intersects the inner surface Si of the vegetative ectoderm T is focused on. The shortest distance between the point P1 and the outer surface So of the vegetative ectoderm T is taken as the thickness of the vegetative ectoderm T in the direction of the radius vector r1. Alternatively, the point P2 where the radius vector r1 intersects the outer surface So of the vegetative ectoderm T is focused on, and the shortest distance between the point P2 and the inner surface Si is taken as the thickness of the vegetative ectoderm T in the direction of the radius vector r1. The calculation of the shortest distance is the same as the calculation process of the thickness of the zona pellucida ZP, and the DistanceTransform function in the OpenCV library can be used.

[0122] Fig.10 4 is a flowchart showing the process of calculating the thickness of the trophectoderm. First, a radial direction is selected (step S401). Specifically, a radial direction is determined by appropriately temporarily setting the values ​​of the deflection angles θ and φ. Next, the position of the inner surface Si of the trophectoderm T in the radial direction in the polar coordinate space ( Fig. 9B The position of the point P1 shown) or the position of the outer surface So ( Fig. 9B The pixel corresponding to the point P2 shown in FIG. 1 is determined as the pixel of interest (step S402). For the following processing, it is preferred to associate and store the distance from the origin O to the pixel of interest with the radial direction determined by the combination of the deflection angles θ and φ. The information representing the distance is represented by the function R (θ, φ).

[0123] Next, the distance between the inner surface Si and the outer surface So of the trophectoderm T in the radial direction is calculated, that is, the thickness of the trophectoderm T (step S403). The thickness can be obtained, for example, by Fig.9A , Fig. 9B Calculate using any of the methods shown.

[0124] The thickness of the vegetative ectoderm T obtained in this way is associated with a radial direction, so it can be expressed as a function of the deflection angles θ and φ. In the following, the thickness of the vegetative ectoderm T in a radial direction determined by the combination of the deflection angles θ and φ is expressed as T(θ, φ). The obtained thickness T(θ, φ) is stored and saved in the memory 37 (step S404). By repeatedly performing the above process for each radial direction (step S405), the thickness distribution representing the thickness of the vegetative ectoderm T in each direction can be obtained.

[0125] exist Figure 3In step S117, a graph showing the spatial distribution of the thickness of the trophectoderm T obtained in this way is produced. As a thickness distribution graph, a two-dimensional graph and a three-dimensional graph can be considered. In addition, as a two-dimensional graph, a method of storing all directions in one graph and a method of dividing the polar coordinate space into several spaces and storing them in each graph can be considered. Specific examples of production of these mapping methods are described below.

[0126] Fig.11A and Fig. 11B is a diagram showing an example of a two-dimensional mapping method. In order to map each direction represented by polar coordinates onto a two-dimensional plane, for example, a quasi-cylindrical projection can be used. Fig.11A As shown in FIG. 1 , as a projection in which the declination angles θ and φ correspond to latitude and longitude, respectively, the same projection as that used to represent the topography of a substantially sphere, that is, the surface of the earth, on a plane map can be used. Several such mapping methods are known, and an appropriate one can be selected according to the purpose.

[0127] Here, a map projection, namely the Sanson projection, is used, which treats the declination angle θ as latitude and the declination angle φ as longitude and can maintain the area of ​​a solid on the map. Pixels having brightness values ​​corresponding to the thickness T(θ, φ) of the trophectoderm T calculated in each radial direction are arranged at each point on the map determined by the combination of the declination angles θ and φ. By indicating the thickness in each direction on the map in this way, the thickness distribution can be visualized. Fig. 11B The thickness of the trophectoderm T at each point is represented by isoluminance lines, instead of directly representing the brightness. It can also be said that the thickness in each direction is represented only by contour lines. In addition, the thickness of the trophectoderm T can also be represented by color or lightness.

[0128] It should be noted that in Fig. 11B In the figure, the shaded area on the upper right indicates the area where the inner cell mass I exists. There is also a vegetative ectoderm T between the inner cell mass I and the zona pellucida ZP, but the shape of the surface in contact with the inner cell mass I does not necessarily represent the shape of each cell. Therefore, in this embodiment, the area where the inner cell mass I exists is excluded from the derivation target of the thickness distribution.

[0129] Another problem in making maps is that when a two-dimensional image represents the surface of a roughly spherical body, all information cannot be accurately represented. For example, in the case of the Samson projection, which is a type of equal-area pseudocylindrical projection, the area of ​​the surface can be accurately represented, but the distance and direction cannot necessarily be accurately represented. In the case of the mapping of this embodiment, the deformation of the peripheral portion of the map plane becomes particularly large.

[0130] Fig. 12A and Fig. 12B is a diagram showing another example of a two-dimensional mapping method. Fig. 12AAs shown in FIG. 1 , the three-dimensional image of the trophectoderm T represented as a substantially spherical surface in the polar coordinate space is divided into four parts by two planes passing through the origin O and being perpendicular to each other. Then, each of the divided subareas Ta to Td is represented on a single figure. At this time, in the figure corresponding to the subarea Ta, the intersection point Ca of the spherical surface and the straight line passing through the origin O and intersecting the two planes at an angle of 45 degrees respectively through the center of the divided spherical surface is substantially the center of the figure. The same is true for the other subareas Tb, Tc, and Td.

[0131] In this way, if Fig. 12B As shown, the trophectoderm T divided into four parts by the polar coordinate space is expressed using only the central part of the drawing plane. Therefore, the drawing has four planes, and the overall thickness distribution of the trophectoderm T can be expressed with high accuracy in each drawing without using the peripheral part where the deformation becomes large.

[0132] It should be noted that, in fact, when considering the division of the sphere corresponding to the trophectoderm T, each partition partially overlaps. Fig. 12B As shown by the middle dotted line, it is more preferable to express the thickness distribution not only in the partition Ta but also in a certain range outside it. The reason is as follows. Although it is divided for the convenience of mapping, the actual vegetative ectoderm T continues to the outside of each partition. By mapping to a certain range outside the partition, the continuity of the thickness distribution of the outer edge of each partition can be expressed.

[0133] In this sense, it is also conceivable to create a plurality of two-dimensional images with different center positions instead of performing segmentation and use only the central portions of these images. The number of segmentations is not limited to the above four but is arbitrary.

[0134] The two-dimensional mapping method described above is suitable for, for example, a bird's-eye view of the overall thickness distribution on a screen or printed paper. However, it is inevitable that the accuracy of the peripheral portion is reduced. Therefore, a three-dimensional map is better in terms of accuracy. In a three-dimensional map, the accuracy of all directions in the polar coordinate space can be guaranteed to be the same.

[0135] When displaying a three-dimensional image on the screen, for example, if a map projected along a certain sight line direction is temporarily displayed and the sight line direction is changed according to the user operation, the convenience of the user can be improved. For example, if an animation display is set in which the sight line direction changes in real time in conjunction with the user operation, the user can perform the evaluation operation as if he or she is directly observing the embryo E in front of him or her.

[0136] As a method of three-dimensional mapping of the vegetative ectoderm T, for example, the following method can be considered. The first method is a method of assigning information indicating thickness distribution in various directions to an approximate spherical surface (or a rotating ellipsoidal surface) indicating the approximate shape of the vegetative ectoderm T. The second method is a method of assigning information on thickness distribution to a curved surface corresponding to the actual three-dimensional shape of the vegetative ectoderm T. In the latter case, it is possible to select a curved surface that corresponds to the inner surface Si of the vegetative ectoderm T or the outer surface So.

[0137] The outer surface So of the vegetative ectoderm T is in contact with the zona pellucida ZP, so its surface shape is relatively smooth. In this sense, the curved surface corresponding to the outer surface So of the vegetative ectoderm T is not substantially different from a spherical surface. On the other hand, the inner surface Si of the vegetative ectoderm T has concavities and convexities corresponding to the shapes of the individual cells constituting the vegetative ectoderm T. Therefore, mapping the curved surface corresponding to the inner surface Si of the vegetative ectoderm T has the advantage of making it easier to imagine the actual concavity and convexity shape.

[0138] It should be noted that as a method of visualizing thickness using a map surface, as in the case of two-dimensional mapping, a method of converting thickness into brightness and expressing it, a method of expressing it by contour lines (isoluminance lines), color distinction, shades, etc., etc., can be used. In the following, examples are given of a case of mapping on a spherical surface and a case of mapping on a curved surface corresponding to the inner surface Si of the trophectoderm T. It is assumed that contour lines are used as a mapping method.

[0139] Fig.13A and Fig. 13B is a diagram showing an example of a three-dimensional mapping method. Fig.13A An example of mapping on a sphere is shown. In addition, Fig. 13B An example of mapping on a curved surface corresponding to the inner surface Si of the vegetative ectoderm T is shown. Here, the three-dimensional diagram is shown in a state of being projected onto a two-dimensional plane, but in fact, mapping on a three-dimensional curved surface is performed. In the case of display output to a screen, it is possible to display in various different sight directions. It should be noted that in these figures, the area indicated by the hatching corresponds to the inner cell mass I.

[0140] For example, Fig.13A The three-dimensional graph can be produced by configuring pixels of brightness corresponding to the thickness calculation result in the direction of the radial vector passing through the position at each position on the spherical surface centered at the origin O in the polar coordinate space. That is, pixels of brightness values ​​corresponding to the thickness T (θ, φ) can be configured at the intersection between the radial vector represented by two deflection angles (θ, φ) and the spherical surface.

[0141] The problem of deformation generated in the two-dimensional diagram is eliminated in the three-dimensional diagram. In the two-dimensional diagram, in addition to the problem of deformation, there is also the problem that the length of the radius vector, that is, the distance from the origin O, is not reflected in the diagram. In particular, when the shape of the sample S deviates greatly from the spherical surface, the concave and convex shown in the diagram may not reflect the actual shape well. In the case of a three-dimensional diagram, such deviation can be reduced by, for example, using a rotating ellipsoid that is closer to the shape of the sample as an approximate curved surface.

[0142] In particular, Fig. 13B As shown in FIG. 1 , by mapping on a curved surface that reflects the actual shape, a more realistic diagram can be produced. The information required for mapping is the position of the inner surface Si (or outer surface So) of the vegetative ectoderm T in each radial direction (the length of the radial vector) and the thickness of the vegetative ectoderm T in that direction. This information is used in the step of calculating the thickness of the vegetative ectoderm T ( Figure 3 In step S116), R(θ, φ) and T(θ, φ) have been obtained respectively. Therefore, a three-dimensional graph is produced by configuring pixels having brightness values ​​corresponding to the thickness T(θ, φ) at positions represented by coordinates (R(θ, φ), θ, φ) in the polar coordinate space.

[0143] It should be noted that Fig. 13B When the three-dimensional graph shown is projected and displayed on a two-dimensional display screen, it is necessary to display the concave-convex shape of the graph surface and the thickness of the trophectoderm T at each position in a recognizable manner. In particular, the shape of the curved surface outside the periphery of the projected two-dimensional graph may be difficult to read from the displayed image. In order to avoid this problem, it is necessary to express the shape of the curved surface and the thickness at each point in a different way. For example, the shape of the curved surface can be expressed by contour lines, and on the other hand, the thickness of each part can be expressed by different shades and colors.

[0144] Back to Figure 3 , the method of counting the number of cells constituting the vegetative ectoderm T based on the thickness distribution mapped in this way is further described. The surface unevenness of the vegetative ectoderm T visualized in the above manner is formed by the shape of each cell. Therefore, in the case where the thickness is represented by brightness, for example, it is considered that the areas with higher brightness than the surroundings scattered on the figure correspond to the protrusions corresponding to each cell in the vegetative ectoderm T. Therefore, if such areas that are thicker than the surroundings, that is, high brightness, are extracted, it is possible to infer where each cell is located. In addition, the number of cells constituting the vegetative ectoderm T can be counted.

[0145] Therefore, the peak appearing in the thickness profile of the trophectoderm T is likely to be a peak indicating a protrusion corresponding to each cell constituting the trophectoderm T. Therefore, in step S118, a clear peak in the thickness profile is searched. In this way, the protrusion is indirectly detected.

[0146] FIG. 14A to FIG. 14C is a diagram showing the principle of peak detection processing. Fig.15 is a flow chart showing the peak detection process. Fig.14A As shown by the solid line in , the relationship between the position in a certain direction and the brightness (i.e., thickness) is plotted along the graph (thickness curve). In this way, there is a change in brightness at each position, so that the brightness stands out relative to the surroundings and a high peak appears. These are candidates for the peak corresponding to the position of the cell.

[0147] First, as described above, an area with a higher brightness than the surrounding area is extracted. Specifically, a maximum value filtering process is performed on the thickness curve shown by the solid line (step S501). By appropriately setting the window size and performing the maximum value filtering process, as shown in FIG. Fig.14A As shown by the dotted line in the middle, the peak width of the profile expands. Through the maximum value processing, the brightness value at the position other than the peak value increases, while the brightness value at the peak value position does not change. Therefore, as described above, the position where the brightness value does not change (indicated by the dotted line in the figure) is extracted as a peak candidate (step S502).

[0148] Fig. 14B Indicates the extracted peak candidates. When the brightness at the peak is extremely high or extremely low, the peak may be a peak different from the detection object, such as image noise, a peak that does not correspond to a cell, etc. Therefore, thresholds Lth1 and Lth2 are set for the brightness value. Peak candidates whose brightness values ​​exceed the threshold Lth1 of the specified upper limit and whose brightness values ​​are lower than the threshold Lth2 of the specified lower limit (step S503) are removed from the peak candidates. Fig.15 This process is described as "upper limit / lower limit process" in Fig. 14B In FIG. 1 , black circles indicate peak candidates that were removed by the process, and white circles indicate peak candidates that were not removed and remained.

[0149] In the purpose of this embodiment to determine each cell constituting the trophectoderm T, the threshold values ​​Lth1 and Lth2 can be set corresponding to the upper limit value and the lower limit value of the appropriate size of the cell constituting the trophectoderm T. In this way, it is possible to exclude extremely large or small peak candidates that cannot be regarded as a cell. As a result, the error caused by the counting result of the number of cells can be reduced. It should be noted that the upper limit value and the lower limit value are set for the peak brightness value here, but either one can be set. In addition, the upper limit / lower limit processing can be omitted as long as it does not cause a large error.

[0150] On the other hand, the cell size is also appropriate in the direction along the surface of the trophectoderm T. In other words, if the peak candidates appropriately represent each cell, the distance between these peaks should be within a specified range corresponding to the size of the cell. Therefore, the distance between mutually adjacent peaks among the remaining peak candidates is calculated (step S504), and in the case of very close peak candidates, one of them is regarded as noise and removed (step S505). In reality, among the peaks having such a relationship, the smaller peak is removed. Fig. 14C An example of the result of removing peak candidates whose distance between adjacent peaks is smaller than a preset minimum distance Dmin is shown.

[0151] In this way, the remaining Fig. 14C The peaks are indicated by white circles in the figure. It is considered that these peak positions represent the positions of the protrusions corresponding to the respective cells. That is, through the processing so far, the protrusions present in the vegetative ectoderm T surrounding the blastocyst cavity B have been detected.

[0152] It should be noted that, in order to make the concept easy to understand, a one-dimensional thickness curve is used to illustrate the principle, but the actual thickness curve is two-dimensional when a two-dimensional graph is used, and is three-dimensional when a three-dimensional graph is used. Therefore, when a two-dimensional graph is used, the maximum value processing in step S501 is two-dimensional filtering. In order to apply a constant window size in all directions, a circular filter is used in this case. In addition, the distance between peaks is also the distance on the plane of the two-dimensional graph.

[0153] In the two-dimensional graph, the information about the length of the radius vector is not reflected, and the deformation becomes larger as it approaches the periphery. Therefore, strictly speaking, it cannot be said that an ideal filtering process is performed in a circular filtering window. However, in practical terms, sufficient accuracy can be obtained for the purpose of detecting the peak value corresponding to each cell position. In particular, Fig. 12A , Fig. 12B As shown in FIG. 1 , in the method of dividing the object into multiple parts and mapping them, the error caused by such a reason can be fully suppressed because the peripheral part of the image where the deformation is large is not used. Fig. 12BThe range shown by the dotted line is at least larger than the size of the filtering window to avoid errors caused by filtering processing.

[0154] In addition, when the thickness curve is represented by a three-dimensional graph, the filter window is a three-dimensional sphere. That is, the maximum value filtering process is a spherical filter. In addition, the distance between peaks is also the distance in the three-dimensional space.

[0155] The thickness of the trophic ectoderm T at each position is projected onto a three-dimensional sphere ( Fig.13A ), the information of the length of the radial vector from the origin O to the pixel of interest is not reflected in the processing. This has the effect of simplifying the processing, but it may also be a cause of error in the filtering process. On the other hand, in the three-dimensional image ( Fig. 13B ), the distance from the origin O is reflected at the position of the pixel of interest. Therefore, by performing spherical filtering processing centered on the pixel of interest, errors can be eliminated and peak detection can be performed appropriately.

[0156] Fig.16A and Fig. 16B is a diagram showing an example of peak positions detected on a two-dimensional graph. Fig.16A As shown by the dotted line, in the region where the thickness of the trophectoderm T is locally prominent compared to the surrounding area, the region satisfying the above conditions is detected as a clear peak (black circle mark). It is believed that these peaks correspond to the cells constituting the trophectoderm T one by one. Fig. 16B This is a graph that plots only the detected obvious peaks on a two-dimensional graph. As described above, when the obvious peaks corresponding to each cell are detected, the number of protrusions is calculated by counting these numbers (step S119). In this way, the number of cells constituting the trophectoderm T can be indirectly known.

[0157] In the region occupied by the inner cell mass I, it is difficult to know the distribution of cells constituting the vegetative ectoderm T. However, by assuming that cells are distributed at the same density as in other regions, the number of cells as the whole of the embryo E can be inferred. That is, by dividing the number of counted cells by the surface area of ​​the vegetative ectoderm T except the part in contact with the inner cell mass I, the cell density is obtained. Then, by multiplying the surface area of ​​the vegetative ectoderm T including the part in contact with the inner cell mass I by the cell density, the total number of cells constituting the whole of the vegetative ectoderm T can be estimated.

[0158] Fig.17A and Fig. 17B: is a diagram showing an example of the peak position detected on the three-dimensional graph. In the three-dimensional graph, as described above, in the region where the thickness of the vegetative ectoderm T locally protrudes compared to the surrounding area, the region that meets the conditions is detected as a clear peak (black circle) corresponding to the cell. By counting it in the entire graph surface, the number of cells constituting the vegetative ectoderm T can be estimated. The processing of the region occupied by the inner cell mass I can also be the same as in the case of the two-dimensional graph.

[0159] As described above, when evaluating cells in a layer covering an internal cavity, quantitative evaluation of cells can be easily performed by detecting protrusions thicker than the surroundings in the thickness curve obtained from the three-dimensional image and treating them as the corresponding parts of each cell. In the case where the outer shape of the object is roughly spherical, it can be easily processed by representing the image data in polar coordinate space.

[0160] It should be noted that the present invention is not limited to the above-mentioned embodiments, and various changes other than the above-mentioned embodiments can be made without departing from the scope of the main purpose. For example, the image processing device 1 of the above-mentioned embodiment has the function of performing OCT photography and optical microscope photography on the sample S and the function of producing and outputting an output image based on the photographing data. However, the image processing method of the present invention can also be executed by a computer device that does not have a photographing function. The computer device only needs to obtain the photographing data obtained by photographing other devices with a photographing function. In order to make this possible, the present invention can also be implemented as a method for enabling a computer device to execute Figure 3 The software program other than step S101 in each processing step.

[0161] The distribution of such a program can be carried out, for example, by downloading it via an electric communication line such as the Internet. Alternatively, it can be carried out by distributing a computer-readable recording medium recording the program. Alternatively, by reading the program into an existing OCT imaging device via an interface, the present invention can also be implemented by the device.

[0162] In addition, for example, in the above-mentioned embodiment, after separating the zona pellucida from the three-dimensional image of embryo E, the remaining structure is further separated into the trophectoderm and the inner cell mass. However, in the gist of the present invention, it is sufficient if the zona pellucida can be separated from other structures, and further separation is not necessary. In addition, the image separated into these regions can also be further processed.

[0163] In addition, for example, in the above-mentioned embodiment, by the three-dimensional image binarization that is obtained by OCT shooting, the structural region of the structure with higher density is separated from the region in addition. By like this, the position of the structure in the image and the occupied range thereof are clarified, thus, for example, there is the advantage of being able to clearly determine the surface position of the embryo when not being affected by the fine structure on the surface, image noise etc. However, it is also possible to perform the separation of the zona pellucida without carrying out this binarization.

[0164] In addition, the above embodiment distinguishes the area occupied by the transparent zone from other areas based on the three-dimensional image obtained by OCT photography, and prompts the user with information useful for evaluation. The information prompt to the user is not limited to this, and it can also be an embodiment that combines various quantitative information obtained by analyzing the three-dimensional image and prompts the user.

[0165] In addition, in the above-mentioned embodiment, various processing is performed after the three-dimensional image data is converted into polar coordinate data. However, even if such coordinate conversion is not performed, the protrusion part can be detected in the same way as described above. For example, when the Distance Transform function is used in the thickness detection, it is not necessary to set the center (center of gravity) of the embryo. Therefore, the thickness can be directly calculated based on the three-dimensional image data of the vegetative ectoderm represented by the orthogonal coordinate system. Then, when using Fig. 13B In the case of the three-dimensional graph shown, the graph can be created in the XYZ orthogonal coordinate space by simply determining the position of the inner surface or the outer surface. Therefore, if these are combined, it can be said that conversion to the polar coordinate system is unnecessary.

[0166] In addition, for example, in the above-mentioned embodiment, the thickness distribution in the vegetative ectoderm is obtained after the zona pellucida and the vegetative ectoderm are separated from each other based on the three-dimensional image of the embryo E. However, if the thickness of the zona pellucida is considered to be constant, even if the thickness distribution is obtained without separating the zona pellucida and the vegetative ectoderm, the thickness variation that appears therein substantially reflects the thickness of the vegetative ectoderm.

[0167] In addition, the above-mentioned embodiment applies the present invention to the evaluation of embryos in the blastocyst stage. However, the scope of application of the present invention is not limited thereto, and can be applied to the observation / evaluation of various cell blocks with cavities inside. For example, the organoid of intestinal epithelial cells has the same cavities as above inside, and its surface is covered by a layer formed by cells. The present invention can also be preferably applied to the evaluation of such cell layers.

[0168] As described above with reference to specific embodiments as examples, in the image processing method of the present invention, for example, the thickness of a cell mass along a radial vector with a point inside the cell mass as the origin can be obtained. When the cell mass has an outer shape similar to a sphere or a rotational ellipsoid, processing can be easily performed by obtaining the thickness in the radial vector direction as described above.

[0169] In addition, for example, based on the three-dimensional image data, polar coordinate data is generated in which the position of each point of the cell mass is represented as a position in the polar coordinate space, and the thickness distribution of the cell mass can be obtained based on the polar coordinate data. In particular, it is preferred to use the center of gravity of the cell mass as the origin of the polar coordinate space. According to this structure, it is easy to process data of a cell mass that is approximately spherical or spheroidal.

[0170] In this case, for example, a cell mass thickness curve can be created that associates the radial direction in the polar coordinate space with the thickness of the cell mass obtained in the radial direction, and based on this, the protrusions can be detected. According to this structure, based on the change in the thickness of the cell mass represented by the thickness curve, it is possible to easily detect the protrusions that protrude relative to the surrounding thickness.

[0171] Specifically, for example, it is possible to perform maximum value filtering on the thickness curve and take the position where the value does not change before and after as the position of the protrusion. Since the thickness of the cell mass at the protrusion is greater than the surrounding area, the thickness value does not change even after the maximum value filtering. On the other hand, in the area with a portion thicker than the surrounding area, the value after the maximum value filtering changes due to the influence of its thickness. By using this, it is possible to detect a portion thicker than the surrounding area.

[0172] For example, a thickness curve can be represented by a three-dimensional graph in which pixels having brightness values ​​corresponding to the thickness of the cell mass at the position corresponding to the inner surface or the outer surface of the cell mass are arranged at the position corresponding to the inner surface or the outer surface of the cell mass. In addition, for example, a thickness curve can be represented by a three-dimensional graph in which pixels having brightness values ​​corresponding to the thickness of the cell mass in the radial direction passing through the position are arranged at each position on the surface of an approximate sphere obtained by approximating the cell mass as a sphere in the polar coordinate space. In these cases, a maximum value filtering process can be performed as a spatial filtering process performed on a spherical region centered on the pixel of interest.

[0173] On the other hand, it is also possible to represent the thickness curve by configuring a two-dimensional graph of pixels having brightness values ​​corresponding to the thickness of the cell mass in one radial direction determined in the polar coordinate space by the two deflection angles represented by the polar coordinate data as coordinate axes on a coordinate plane. In this case, as a spatial filtering process for a circular area centered on the pixel of interest on the two-dimensional graph, a maximum value filtering process can be performed. It should be noted that as a method of representing the thickness distribution in three-dimensional space by a two-dimensional graph, for example, a quasi-cylindrical projection can be used.

[0174] In addition, for example, the polar coordinate space can be divided into multiple planes passing through the origin, and the two-dimensional map and the protrusion detection can be performed separately for each of the divided spaces. When the shape of an object in a three-dimensional space is represented by a two-dimensional map, deformation of the object on the map cannot be avoided. By dividing the space and making separate maps, large deformations can be avoided at certain locations, thereby improving the detection accuracy of the thickness distribution.

[0175] In addition, for example, the image processing method of the present invention may also include a step of measuring the number of protrusions detected. By obtaining quantitative information in this way, more useful information can be provided to the user. For example, when the cell mass of the object is an embryo in the blastocyst stage or an intestinal epithelial cell organoid, the present invention can be applied to count the cells constituting its surface.

[0176] In addition, in an image processing method that uses an image of a fertilized egg as a processing object, for example, a classification algorithm that has been pre-machined using an image of a zona pellucida taken with an optical microscope as a teacher image can be used to determine an area corresponding to the zona pellucida. According to such a structure, by using an appropriate classification algorithm, an area having a high morphological feature of the zona pellucida can be determined with high accuracy from an image taken with an optical microscope.

[0177] For example, as a classification algorithm, a semantic segmentation method can be used. According to this method, an image can be segmented using pixel units based on its features. Therefore, the area of ​​the zona pellucida can be determined with high precision from a microscope image, and its thickness can be accurately evaluated.

[0178] In addition, for example, the optical microscope image of the fertilized egg may be taken multiple times with different focal depths, and the average thickness of the zona pellucida may be obtained from the image with the largest area of ​​the region corresponding to the zona pellucida among these images. The image that determines the zona pellucida in a wide range is most likely to be the image that is most focused on the zona pellucida, that is, the image that is taken in a state where the zona pellucida is clear. By using such an image, the thickness of the zona pellucida can be obtained with high accuracy.

[0179] In this case, for example, the area of ​​the region determined using the classification algorithm can be multiplied by a coefficient determined based on the magnitude of the brightness change on the boundary between the inside and outside of the region, as a value representing the area of ​​the region corresponding to the transparent band. For example, the square of the difference in brightness between the inside and outside of the boundary can be used as the coefficient. With such a configuration, the risk of misjudging an image whose apparent transparent band area becomes larger due to shaking during shooting or the like as a focused image can be reduced.

[0180] In addition, for example, the region in the structural region that is within the range of average thickness from the outer edge can be regarded as the zona pellucida. In a properly cultured fertilized egg, the thickness of the zona pellucida covering its surface is approximately uniform. Therefore, when observing from the outer edge of the fertilized egg in the three-dimensional image, it can be said that the region up to the depth corresponding to the average thickness obtained from the microscope image is the region occupied by the zona pellucida.

[0181] In this case, the region outside the zona pellucida in the structural region may be further regarded as the region of the inner cell mass and the vegetative ectoderm. The main structural bodies constituting the fertilized egg in the blastocyst stage properly cultured are the zona pellucida, the vegetative ectoderm and the inner cell mass. Therefore, it can be considered that the vegetative ectoderm and the inner cell mass occupy the region outside the zona pellucida in the structural region.

[0182] In addition, for example, an area occupied by pixels having a brightness value greater than a specified value may be determined based on the three-dimensional image data, and the area may be regarded as a structural area. In an OCT image, a higher-density structure has a higher brightness, and the surrounding low-density area has a lower brightness. Therefore, if an area having a certain degree of brightness is determined, it is highly likely that the area is a structural area occupied by a high-density structure.

[0183] Industrial Applicability

[0184] The present invention is suitable for assisting the observation and evaluation of cell masses having cavities inside. For example, it can be applied to assisting the evaluation of the state of cultured embryos (fertilized eggs) and selecting good embryos with a higher pregnancy success rate in reproductive assisted medicine.

Claims

1. An image processing method, wherein an image of a cell mass having a cavity inside is taken as a processing object, in, include: a step of performing optical coherence tomography on the cell mass and acquiring three-dimensional image data representing a three-dimensional image of the cell mass; A step of acquiring a thickness distribution indicating the thickness of the cell block between an outer surface facing the external space and an inner surface facing the cavity in the cell block at each position based on the three-dimensional image data; and a step of detecting a protruding portion having a thickness greater than that of a surrounding area based on the obtained result of the thickness distribution, Based on the three-dimensional image data, polar coordinate data is generated, which represents the position of each point of the cell mass as a position in a polar coordinate space; Based on the polar coordinate data, the thickness distribution of the cell mass is obtained. detecting the protrusion part based on a thickness curve of the cell mass that associates the radial direction of the polar coordinate space with the thickness of the cell mass obtained in the radial direction, A maximum value filtering process is performed on the thickness curve, and a position where the value does not change before and after the filtering process is used as the position of the protrusion.

2. The image processing method according to claim 1, wherein: The thickness of the cell mass is determined along a radial vector having a point inside the cell mass as an origin.

3. The image processing method according to claim 1 or 2, wherein: The center of gravity of the cell mass is the origin of the polar coordinate space.

4. The image processing method according to claim 1, wherein: representing the thickness curve by a three-dimensional graph produced by arranging pixels having brightness values ​​corresponding to the thickness of the cell mass at a position corresponding to the inner surface or the outer surface, The maximum value filtering process is performed as a spatial filtering process on a spherical area centered on the pixel of interest.

5. The image processing method according to claim 1, wherein: The thickness curve is represented by a three-dimensional graph, the three-dimensional graph being produced by arranging pixels having brightness values ​​corresponding to the thickness of the cell mass in a radial direction passing through the position at each position on the surface of an approximate sphere obtained by spherically approximating the cell mass in the polar coordinate space, The maximum value filtering process is performed as a spatial filtering process on a spherical area centered on the pixel of interest.

6. The image processing method according to claim 1, wherein: representing the thickness curve by a two-dimensional graph produced by arranging pixels having brightness values ​​corresponding to the thickness of the cell mass in a radial direction determined by the two deflections represented by the polar coordinate data in a coordinate plane having two deflections as coordinate axes in the polar coordinate space, The maximum value filtering process is executed as a two-dimensional filtering process on a circular area centered on the pixel of interest on the two-dimensional map.

7. The image processing method according to claim 6, wherein: The two-dimensional image is produced by quasi-cylindrical projection.

8. The image processing method according to claim 6, wherein: The polar coordinate space is divided into a plurality of parts by a plane passing through the origin, and the creation of the two-dimensional map and the detection of the protrusions are performed separately for each of the divided spaces.

9. The image processing method according to claim 1 or 2, wherein: The image processing method further includes a step of counting the number of the detected protrusions.

10. A computer program product, comprising a computer program, wherein when the computer program is executed by a processor, the computer program implements the following steps: A step of performing optical coherence tomography on a cell block having a cavity therein, and acquiring a thickness distribution indicating the thickness of the cell block between an outer surface facing an external space and an inner surface facing the cavity at each position based on three-dimensional image data indicating a three-dimensional image of the cell block; and a step of detecting a protruding portion having a thickness greater than that of a surrounding area based on the obtained result of the thickness distribution, Based on the three-dimensional image data, polar coordinate data is generated, which represents the position of each point of the cell mass as a position in a polar coordinate space; Based on the polar coordinate data, the thickness distribution of the cell mass is obtained. detecting the protrusion part based on a thickness curve of the cell mass that associates the radial direction of the polar coordinate space with the thickness of the cell mass obtained in the radial direction, A maximum value filtering process is performed on the thickness curve, and a position where the value does not change before and after the filtering process is used as the position of the protrusion.

11. A computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of the image processing method according to claim 1.

12. An image processing method, in, include: A step of acquiring three-dimensional image data obtained by optical coherence tomography and two-dimensional image data obtained by optical microscopy of a fertilized egg at the blastocyst stage; A step of determining a structural region corresponding to the structural body of the fertilized egg based on the three-dimensional image data; A step of determining a region corresponding to the zona pellucida of the fertilized egg based on the two-dimensional image data and calculating an average thickness thereof; and The structure region is divided into a region within the range of the average thickness from the outer edge toward the inside of the structure region and other regions.

13. The image processing method according to claim 12, wherein: The area corresponding to the transparent zone is determined by using a classification algorithm that performs machine learning in advance using an image of the transparent zone taken with an optical microscope as a teacher image.

14. The image processing method according to claim 13, wherein: The semantic segmentation method is used as the classification algorithm.

15. The image processing method according to claim 13 or 14, wherein: The optical microscope image of the fertilized egg is taken a plurality of times while changing the focal depth, and the average thickness of the zona pellucida is determined from the image having the largest area of ​​the region corresponding to the zona pellucida among the images.

16. The image processing method according to claim 15, wherein: The area of ​​the region determined using the classification algorithm is multiplied by a coefficient determined based on the magnitude of brightness change at the boundary between the inside and outside of the region, as a value representing the area of ​​the region corresponding to the transparent zone.

17. The image processing method according to claim 16, wherein: The coefficient is the square of the difference between the brightness values ​​inside and outside the boundary.

18. The image processing method according to any one of claims 12 to 14, wherein: The region in the structural region that is within the range of the average thickness from the outer edge is regarded as the transparent zone.

19. The image processing method according to claim 18, wherein: The region outside the zona pellucida in the structural region is regarded as the region of the inner cell mass and vegetative ectoderm.

20. The image processing method according to any one of claims 12 to 14, wherein: Based on the three-dimensional image data, a region occupied by pixels having a brightness value greater than a predetermined value is determined, and the region is regarded as the structure region.

21. A computer program product, comprising a computer program, wherein when the computer program is executed by a processor, the computer program implements the following steps: The step of determining a structural region corresponding to the structural body of the fertilized egg based on three-dimensional image data obtained by optical coherence tomography of the fertilized egg at the blastocyst stage; A step of determining a region corresponding to the zona pellucida of the fertilized egg based on two-dimensional image data obtained by photographing the fertilized egg under an optical microscope and calculating an average thickness thereof; and The step of dividing the structure region into a region within the range of the average thickness from the outer edge toward the inside of the structure region and other regions.

22. A computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of the image processing method according to claim 12.

Images