System and method for subpopulation identification

Abstract

The present invention provides systems and methods for identifying and distinguishing subpopulations within a particle mixture. The systems and methods enable continuous calibration of particle classification in a mixture to provide consistency of operation and reduce inter-batch and intra-batch processing variations. The systems and methods produce advantageous sorted particle products.

Description

Background Art

[0001] Often this application is related to flow cytometry. Specifically, this application involves systems and techniques to determine whether a particular cellular activity is included in or excluded from a subpopulation or a gate.

[0002] Flow cytometers can use light detection to evaluate the characteristics of particles (such as cells) passing through the flow cytometer. In some applications, the cytometer can detect light emitted by cells, including fluorescence, light emitted by DNA intercalating dyes. The ability to evaluate and classify particles is fundamental to the operation of the cytometer.

[0003] Cytometers measure characteristics of particles and provide information about the measured characteristics. Typically, this information is presented in a graphical form such as a histogram or scatter plot, although it is also commonly collected and stored in a computer-readable form. Information related to these measurements can be used to identify subpopulations within the initial collection of particles, or to define a subjective or objective value for the characteristic being measured. This is typically done by defining a gate using the graphical representation. Generally, the process of drawing gates around populations is done manually: an operator visually inspects a graphical representation of a mixture of particles as the mixture is analyzed and manually sets the points used to distinguish between populations or subsets of the mixture. Once subsets have been defined by the drawing of gates, all subsequent events can be defined as falling inside the gate (i.e., belonging to the defined subset) or outside the gate (i.e., not belonging to the defined subset and / or belonging to another subset). This has important implications for many applications, including quantitative flow cytometry and cell identification such as fluorescence-activated cell sorting (FACS).

[0004] One application is to determine whether a sperm cell has an X chromosome (which can produce a female fertilized egg) or a Y chromosome (which can produce a male fertilized egg). A sperm with two X chromosomes may have about 3% more DNA than a sperm with one X chromosome and one Y chromosome. By identifying the chromosomal content of sperm cells, it is possible to produce a relatively high purity population of sperm cells with either X chromosomes or Y chromosomes. For example, such a population can be produced by retaining the desired sex and killing the undesirable sex (or by isolating the two sexes). Large, highly pure populations may be used to inseminate female animals (such as cattle, horses, sheep, goats, pigs, dogs, cats, camels, elephants, bulls, buffaloes, or the like) to obtain the desired male or female offspring with a relatively high probability.

[0005] The characteristics of the chromosomes can be identified by staining the DNA of the sperm cells with a fluorescent dye (e.g., a DNA intercalating dye). The stained sperm cells may be forced to flow in a narrow stream or band and passed through an excitation or irradiation source, such as a laser beam. When the stained sperm cells (multiple particles) are irradiated, the fluorescein in the multiple particles will emit responsive fluorescence. The number of fluorescent lights can be varied at least in part based on the relative number of at least one particle differentiation feature in each of the multiple particles (e.g., the relative number of DNA indicating the presence of an X or Y chromosome). The fluorescence can be received by one or more optical elements, which ultimately focus the received light onto a detection component, such as a photomultiplier tube ("PMT") or an avalanche photodiode (APD). The detection component can generate an electrical analog signal in response to the received light. The analog signal can vary depending on the amount of light received. The signal can then be processed (e.g., digitized and analyzed by a processor) to assess the chromosomal content of the sperm cells, which can be graphically represented as a histogram (e.g., representing the data as a vector of peak voltage values, indexed as a function of bin width) or a scatter plot (e.g., mapping peak voltage values ​​and regional voltage values ​​into a vector for displaying the most recent data in the form of a spectrogram). These graphical representations can provide a basis for defining subpopulations and / or drawing gates. This is typically accomplished by the user using her best judgment to select points for establishing gate or threshold boundaries.

[0006] Current systems for setting limits or thresholds to define subpopulations suffer from several significant deficiencies. First, the selection of defining criteria (i.e., setting the values ​​of the gates or threshold points) is performed by individual operators. This is fundamentally subjective and can be variable between operators, and even between samples run by the same operator, or between samples run by the same operator. This ultimately leads to variability in the operation of the flow cytometer system. In the case of sex sorting cytometers, this can result in significant differences in the yield (i.e., the number of viable cells collected) or purity (i.e., the proportion of cells collected that have an X or Y chromosome) of the sexed sperm cell product.

[0007] Furthermore, current subpopulation definitions and differentiation methods result in a differentiation benchmark that is essentially static during operation—the gate or threshold remains constant during operation. In practice, this means that an operator must monitor the operation of the subpopulation during operation, and / or the relationship of the gate or threshold to the subpopulation may change during operation due to a variety of factors, such as changes in temperature or pressure, or loss of signal.

[0008] The systems and techniques described herein overcome these and other disadvantages and provide significant benefits that will be appreciated by those skilled in the art.

[0009] summary

[0010] According to an embodiment of the invention, the present invention provides a system for identifying subpopulations in a particle mixture. According to a specific embodiment, the system includes a sample delivery device, a cytometer, a control device and a sorting device. The control is configured to implement a self-calibration recognition or identification technology. The self-calibration recognition or identification technology uses information obtained from particles in the mixture to define subpopulations in the mixture. The system then classifies the particles in the mixture into subpopulations on this basis. The definition of the subpopulations is continuously refined by integrating additional information obtained from additional particles in the mixture. The above-mentioned control equipment is also configured to receive event data from the cytometer in the form of values ​​corresponding to the particles detected in the sample, distinguish the measured characteristic values ​​of each subpopulation in the sample, assign a critical value range to each sample based on the distinguished subpopulations, update the critical range of values ​​based on the measured characteristic values ​​that distinguish the subpopulations in the second group of particles, and classify the particles in the mixture into subpopulations based on the updated critical value range.

[0011] According to certain embodiments of the invention, the control can be further configured to identify at least one critical value for a first feature from a distribution; obtain a value of a second feature related to the critical value of the first feature; generate a distribution of values ​​of the second feature; identify at least one critical value of the second feature from the distribution of values ​​of the second feature; co-locate the critical values ​​of the first feature and the second feature in a feature space; and determine at least one discriminant value from the critical values ​​co-located in the feature space.

[0012] According to certain inventive techniques, a method includes: identifying subpopulations within a particle mixture, including: measuring characteristics of a first group of particles in a particle mixture having two or more subpopulations, distinguishing each subpopulation by the characteristics; distinguishing the measured characteristic values ​​of each subpopulation in the mixture; assigning a critical value range to each subpopulation based on the distinction; after assigning the critical value range, measuring characteristics of a second group of particles in the particle mixture; after measuring the characteristics of the second group of particles in the particle mixture, distinguishing the measured characteristic values ​​of each subpopulation in the second group of particles; updating the critical range of values ​​based on the distinction of the measured characteristic values ​​of each subpopulation in the second group of particles; and assigning subpopulation classifications to the particles in the mixture based on the updated critical value range.

[0013] According to certain inventive techniques, the method further comprises classifying the particle according to the classification assigned to the particle based on the updated threshold range. According to certain inventive techniques, the method further comprises destroying or losing the ability of the subpopulation assigned to the classification. According to certain inventive techniques, the particle is a sperm cell. According to certain inventive techniques, the sperm cell is a bovine sperm cell. According to certain inventive techniques, the sperm cell is characterized by having an X chromosome or a Y chromosome. According to certain inventive techniques, the sperm cell is labeled with a fluorescent DNA binding dye.

[0014] According to certain inventive techniques, a method includes identifying multiple groups of particles in a mixture of particle groups, each subgroup having a characteristic different from other subgroups, including: providing a mixture of particles; detecting a characteristic of a first group of multiple particles in the mixture to generate a signal indicating the characteristic of the first group of particles; generating event data for the signal, the signal indicating the characteristic of the first group of particles; compiling the event data of the first group of particles; generating a distribution of the event data for the first group of particles; calculating at least one discrimination value defining a subgroup of particles based on the distribution of the event data of the first group of particles; and classifying the particles based on the discrimination value defining the subgroup. According to certain inventive techniques, the method also includes including or excluding events based on the classification indication classification device. According to certain inventive techniques, the method also includes destroying or incapacitating the subgroup assigned the classification. According to certain inventive techniques, the particles are sperm cells. According to certain inventive techniques, the sperm cells are bovine sperm cells. According to certain inventive techniques, the sperm cells are characterized by having an X chromosome or a Y chromosome. According to certain inventive techniques, the sperm cells are labeled with a fluorescent DNA binding dye.

[0015] According to certain inventive techniques, calculating at least one discriminant value includes at least defining a critical value associated with the detected feature. According to certain inventive techniques, calculating at least one discriminant value also includes: identifying at least one critical value for a first feature from a distribution of event data for a first group of particles; obtaining a value of a second feature associated with the critical value of the first feature; generating a distribution of values ​​of the second feature; identifying at least one critical value for the second feature from the distribution of values ​​of the second feature; co-locating the critical values ​​of the first and second features in a feature space; determining at least one discriminant value from the critical values ​​co-located in the feature space. According to certain inventive techniques, calculating at least one discriminant value also includes adjusting the critical value to include event data from a second group of particles. According to certain inventive techniques, the critical value includes one or more peak characteristic values ​​for a first subpopulation, and one or more peak characteristic values ​​for a second subpopulation. According to certain inventive techniques, the peak characteristic values ​​include peak fluorescence intensity and peak total fluorescence. According to certain inventive techniques, the critical value also includes a valley value between the peak characteristic values ​​of the first and second subpopulations. According to certain inventive techniques, the discriminant value includes a gate.

[0016] According to certain inventive techniques, a sexed sperm cell sample is produced by any of the methods of any of the preceding claims.

[0017] The foregoing summary and the following detailed description of certain technologies and systems of the present application will be better understood when read in conjunction with the accompanying drawings. For illustrative purposes, certain systems and technologies are shown in the accompanying drawings. However, it should be understood that the claims are not limited to the arrangements and instrumentalities shown in the accompanying drawings. Summary of the invention

[0018] Certain inventive systems and techniques discussed herein relate to a cytometer that includes a mechanism that uses known relationships in characteristic frequencies or quantum relationships of particles in a mixture to identify, define, or distinguish subpopulations in a mixture. Advantageous systems and techniques optimize throughput (yield and / or purity), reduce variability and increase operational consistency, and reduce manual user input and variability associated therewith. These advantageous systems and techniques provide decision schemes that reliably distinguish particles having different characteristics, thereby providing different classifications.

[0019] The systems and methods of the present invention can be used to identify one or more subpopulations of particles in a mixture. This can be achieved by detecting or evaluating a single feature or multiple features of the particles. Features can be evaluated using different identification means, including, for example, fluorescent markers (e.g., fluorescent dyes or fluorescent markers coupled to one or more target-specific molecules, such as antibodies or antibody fragments or other molecular markers). On the other hand, the methods and systems of the present invention include sorting cells, such as sperm cells, by distinguishing differences in DNA content. The cells are stained with DNA-specific dyes by mixing sperm cells with dyes and then electroporating. Sperm cells can be maintained at a temperature that increases sperm cell motility, typically at or below 39°C. Sperm cells are then passed through an excitation light source to cause the stained DNA to emit fluorescence, and then through a tool for detecting fluorescence and fluorescence and a device for cell sorting, wherein cells are sorted by DNA content, and sorted sperm cells are collected. The method and apparatus are suitable for sorting mammalian sperm cells, such as sperm cells from cattle, pigs, rabbits, alpacas, horses, dogs, cats, ferrets, rats, mice, and buffaloes. Both membrane-permeable dyes and non-membrane-permeable dyes can be used. Useful dyes include dyes from the SYTOX blue, orange and green series, cyanine dimers and monomers, POPO-1, BOBO-1, YOYO-1, TOTO-1, JOJO-1, POPO-3, LOLO-1, BOBO-3, YOYO-3, TOTO-3, PO-PRO-1, BO-PRO-1, YO-PRO-1, TO-PRO-1, JO-PRO-1, PO-PRO-3, LO-PRO-1, BO-PRO-3, YO-PRO-3, TO-PRO-3, TO-PRO-5, acridine homodimers, 7-aminoactinomycin D, ethidium bromide, homodimeric ethidium-1, homodimeric ethidium-2, ethidium nonazide, nuclear yellow, propidium iodide. Other useful dyes include blue, green, orange and red fluorescent dyes from SYTO 40, Hoechst dyes and dihydroethidium. To enhance the signal, nanoparticles such as quantum dots and metal nanoparticles can be introduced. These particles can be labeled with targeting molecules. The sorting efficiency can be greater than 90%, while the sperm cell survival rate can be greater than 30%, typically greater than 90%. In addition, the dye is introduced into the cells through an osmotic gradient rather than electroporating the cells. The cells are first incubated under hypertonic conditions and then transferred to hypotonic conditions; DNA staining dyes can be added to one or both of the hypertonic and hypotonic solutions. After the cells are stained, they can be sorted or further processed. In other aspects, the systems and methods of the present invention may involve detecting differences in DNA content without the use of DNA stains, for example, U.S. Patent Nos. 8,502,148, 9,835,552; and 9,683,922.On the other hand, the systems and methods of the present invention may involve the use of other analytical means to detect characteristics of the particles, including, for example, side scatter, forward scatter or back scatter of light, or live / dead stains (such as propidium iodide or other active stains).

[0020] In another aspect, the systems and methods of the present invention may involve establishing or modifying a threshold or gate. A threshold generally refers to a value, usually associated with a specific feature, that can be used to classify a single event. For example, an intensity threshold can be set to a value that allows for the distinction of fluorescence that indicates a live cell (e.g., above a threshold) or a dead cell (below a threshold). Similarly, a gate refers to a set of points or vertices that define an area in a two-dimensional graphical representation of a detected feature. For example, a gate may refer to a set of four or more points that are connected to form a polygon around a subpopulation. A gate can be used to define or distinguish a specific subpopulation. The individual points or vertices that make up the gate can be discriminant values ​​calculated from the event data. In a scatter plot, each point or vertex of the gate has X and Y coordinates. The values ​​of these coordinates can be determined by calculation based on the compiled event data and / or by the values ​​contained in the parameter set. The gate can then be used to classify individual events to distinguish whether the event belongs inside or outside the gate. For example, a gate can be determined to include a specific cell type based on the presence of a fluorescent marker that binds to two cell surface markers; cells with both markers will fall within the gate, while cells lacking one or both markers will fall outside the gate.

[0021] Sorting system

[0022] On the one hand, a system for optimizing particle identification and classification is provided. On the one hand, the system includes a sample delivery device. A mixture of particles is delivered to a cytometer system by a sample delivery system. As will be appreciated by those skilled in the art, sample delivery can be accomplished using a syringe pump, a peristaltic pump, and / or increasing the pressure difference between the sample and the system. Typically, the sample delivery system is capable of discretely or continuously sampling a mixture of particles in a liquid suspension and moving the mixture of particles to the cytometer. In some applications, the sample delivery system may also include mixing the mixture of particles with another fluid, and / or delivering a buffer or sheath fluid. The sperm sample to be sorted may be a freshly collected sample from an animal source (e.g., a cow, horse, pig, or other mammalian source), or a thawed, previously frozen sample. In addition, the sample may be a single semen from the same mammal, multiple semens, or multiple pooled semens from two or more animals.

[0023] Various collection methods are known, including the gloved method, the use of an artificial vagina, and electro-ejaculation. The sperm is collected perfectly or quickly transferred to an insulated container to avoid rapid temperature changes from physiological temperature (usually about 35°C to about 39°C). Depending on the species and the particular animal, each milliliter of semen typically contains about 0.5 to 15 billion sperm.

[0024] Regardless of the collection method used, aliquots can be extracted from the sperm sample and evaluated for various characteristics, such as sperm concentration, sperm motility, sperm progressive motility, sample pH, sperm membrane integrity, and sperm morphology. These data can be obtained by sperm examination according to standard and well-known procedures, for example, using a Hamilton-Thorne motility analyzer (IVOS) (see, Farrell et al. Theriogenology (1998) 49 (4): 871-9; and US Patent Nos. 4,896,966 and 4,896,967).

[0025] In another aspect, the system of the present invention includes a cytometer. The cytometer detects features of particles as events and provides signals corresponding to the events to a processor. According to certain aspects, a flow cytometer device includes: a flow chamber configured to guide a fluid flow including sample particles through a particle interrogation position; and a laser configured to emit electromagnetic radiation along a beam path to the particle interrogation position; a detector configured to receive electromagnetic radiation from the interrogation position; output a time-varying analog signal indicating the intensity of the received electromagnetic radiation; at least one amplifier configured to amplify the time-varying analog signal; an analog-to-digital converter configured to receive the amplified time-varying analog signal and generate a corresponding digitized output signal. Further, the detector can be a photomultiplier tube (PMT) or an avalanche photodiode (APD). The cytometer further includes a temperature sensor configured to sense temperature and generate a corresponding temperature signal encoding temperature data, wherein the voltage adjustment circuit is further configured to adjust the voltage of the feedback loop based at least in part on the temperature data and the voltage measured between the filter and the avalanche photodiode. The voltage regulation circuit may be further configured to regulate the voltage on the feedback loop based at least in part on at least one characteristic of the avalanche photodiode, temperature data, and a voltage measured between the filter and the avalanche photodiode. The voltage regulation circuit may be further configured to regulate the voltage on the feedback loop based at least in part on at least one temperature data characteristic of the avalanche photodiode and a voltage measured between the filter and the avalanche photodiode. The voltage regulation circuit may also include: an analog-to-digital converter configured to convert the voltage measured between the filter and the avalanche photodiode into a digital measurement signal encoding the measured voltage data; and a processor configured to at least process the measured voltage data to generate a digital adjustment signal; a digital-to-analog converter configured to convert the digital adjustment signal into an adjustment voltage, wherein the adjustment voltage affects the voltage on the feedback loop. The processor may be configured to process at least the temperature data and the measured voltage data to generate a digital adjustment signal. The processor may be configured to process at least one temperature data, data corresponding to at least one characteristic of the avalanche photodiode, and the measured voltage data to generate a digital adjustment signal. The at least one characteristic of the avalanche photodiode may include at least one of a breakdown voltage and a reverse bias voltage corresponding to a predetermined optical gain. The cytometer may further include: a first amplifier configured to amplify a voltage at an anode of the avalanche photodiode to form a first amplified voltage; and a second amplifier configured to amplify the first amplified voltage to generate a second amplified voltage. The power supply may include a DC / DC power supply.The avalanche photodiode may be: arranged to receive an amount of fluorescence emitted by each of a plurality of particles; the received amount of fluorescence may vary based at least in part on the relative amount of at least one particle differentiation characteristic present in each of the plurality of particles; the avalanche photodiode may be configured to convert the received amount of fluorescence light into at least one signal that varies based on the received amount of fluorescence light.

[0026] A stream of particles can be passed through a flow chamber so that they can form a single file at certain key locations. The flow chamber can guide a fluid including a mixture of particles to flow through a particle interrogation position. The mixture of particles may have been previously dyed with a dye, such as a DNA intercalating dye. When this dye produces a response light in response to exposure to a light source (or electromagnetic radiation), it can emit fluorescence (or cause fluorescence). As the particles pass through one by one, they can be exposed to a beam of electromagnetic radiation (e.g., light of a given wavelength) generated by a laser and emitted to the particle interrogation position along a beam path (and optionally associated optical devices). Such associated optical devices can include lenses, filters, etc. Particles (including particles dyed with dyes) exposed to electromagnetic radiation can cause response events, including light scattering, absorption and / or fluorescence emitted by the dyed particles. In some aspects, the amount of light generated by the fluorescence can vary according to the detectable degree, depending on whether the sperm particle carries an X chromosome or a Y chromosome.

[0027] The light generated by the fluorescence can be received by an optical device (e.g., a lens, a filter, etc.) and can be focused onto the active area of ​​the photodetector. According to certain inventive techniques, the photodetector may include an APD or a PMT. The photodetector may generate an output signal corresponding to (linearly or nonlinearly varying) the amount of electromagnetic radiation (e.g., light of a given wavelength) received from the interrogation location. The photodetector output signal may include a time-varying analog signal indicating the intensity of the received electromagnetic radiation. The photodetector output signal may ultimately be transmitted to a processor (which may include one or more processors that control part or all of the operation of a flow cytometer or sorting device). The photodetector output signal may be amplified by at least one amplifier (or two or more amplifiers in series) before being transmitted to the processor. In addition, the photodetector output signal (e.g., an amplified signal) may be digitized before being transmitted to the processor.

[0028] In another aspect, the system may include a processor that receives a signal and generates event data. The processor may further compile the event data into one or more data sets. The processor may also distribute the event data in the data sets, perform discriminant value calculations; classify events according to the discriminant values ​​and make sorting decisions; and / or provide instructions to a sorting device. The processor may include a single component or may include multiple components. The processor is configured to analyze the digitized output signal received from the flow cytometer. Generally speaking, the electrical signal from the cytometer detector (i.e., PMT or APD) is converted into digital information by an A / D converter, and the A / D converter provides the corresponding digital information to the microprocessor. In response to the information, the microprocessor controls the sorting system.

[0029] The electrical signal output from the photodetector of the cytometer system is a time-varying analog voltage signal indicating the amplitude of the fluorescence emitted by each particle at any time under the irradiation of the laser beam. Therefore, the analog signal (also referred to as analog output) is in the shape of a time-varying waveform pulse. Generally, a waveform pulse is defined as a waveform or a part of a waveform containing one or more pulses or a portion of a pulse. Therefore, the amplitude of each waveform pulse at any time represents the relative rate of photon emission of each cell at that time when the cell passes through the laser beam. In one embodiment, X chromosome bovine sperm cells have a higher DNA content (e.g., about 3.8%) than Y chromosome bovine sperm cells. As a result, live X cells labeled with fluorescent dye as described above will produce waveform pulses 497 that are different from any other labeled cells. By analyzing the pulses as described below (see signal processing, slit scanning, and critical slope difference), each cell can be identified as an X cell or not identified as an X cell (~X). In general, as used herein, X cells refer to live X cells, Y cells refer to live Y cells, and ~X cells refer to a combination of live Y cells and cells that otherwise produce detectable fluorescence emission but cannot be determined as live X cells with a reasonable probability.

[0030] The timing of each waveform pulse indicates the location of each cell in the flow. Since the rate at which particles are transported through the flow cytometer remains constant, and since the distance between the detector and the sorting device is known, the location of each particle is known. Therefore, the microprocessor can calculate the instant when the sorting event occurred (i.e., when the ablation laser fired). Since the microprocessor receives event data for each particle, including the characteristics of the particle, the microprocessor can track (or enumerate) the number of particles that have been sorted in each subpopulation. Based on the parameter settings of the control device and the classification strategy (see below), the microprocessor determines how to sort the particle mixture.

[0031] In another aspect, the system includes a control, in the form of a microprocessor (or other digital or analog control and / or processor or combination thereof), controlling the discrimination of particles in a cell mixture. In a further aspect, the control receives information from a cytometer device, including measured characteristic data of particles in a particle mixture. The control also includes information that a parameter set exists. Using these information sources, the control is able to distinguish the measured characteristic values ​​of subpopulations in a group of particles; assign a critical value range to two or more subpopulations based on the distinguishable values ​​of the measured characteristics; measure a second group of particles in the particle mixture and distinguish the values ​​of the measured characteristics of subpopulations in the second group of particles; use the values ​​of the measured characteristics of the subpopulations distinguished for the second group of particles to update the critical range assigned to the subpopulation; and classify the particles in the mixture according to the subpopulation classification using the updated critical value range. The controller provides an output signal to control the sorting system according to the classification. The control can be implemented within a single microprocessor. In addition, some or all functions can be integrated into one or more processors. In addition, the processing can be implemented by analog circuits or a combination of analog and digital circuits.

[0032] The control may provide output signals to other parts of the cytometer system. In addition, the control may be adapted to process information and provide output signals in real time. Broadly speaking, the term "real time" refers to operation of a control device that matches human perception of time or operation of one or more processors of the control that matches the rate of operation of an associated physical or external process. In certain contexts, the term "real time" may mean that the system reacts to events before the event is out of date.

[0033] On the other hand, the system of the present invention includes a sorting device. The sorting device performs a sorting operation on particles in a particle mixture according to instructions received from the control device. The sorting device can be any type of sorting machine known to those skilled in the art, including, for example, a droplet sorter, a mechanical sorter, a laser ablation ("zapper") system, or any other technology known to those skilled in the art. In one embodiment, as described in U.S. Pat. No. 9,588,100, the sorting system includes a focused energy device. The focused energy device is used to provide focused energy pulses to the particles. The focused energy device can be a thermal, electrical, optical or electromagnetic device that will have the required wavelength and will deliver high peak power to the target particles at a very high repetition rate (or pulse frequency). In one embodiment, the focused energy device is triggered at a predetermined time (i.e., milliseconds) after the controller is activated (the time is set according to the speed at which the particles move through the channel) and pulses are emitted to the selected or targeted (i.e., unwanted) particles.

[0034] Pulsed lasers include, for example, mode-locked, Q-switched, and those lasers that use both mode-locked and Q-switched techniques. For example, focused energy devices such as Avia 355-5-100 (manufactured by Coherent, Inc., Santa Clara, Calif.) or Explorer XP lasers Q-switched lasers from SpectraPhysics can operate in pulse-on-demand mode and can provide 15ns or less energy pulses to a target object at a rate of more than 1000 pulses per second. In one embodiment, the pulse energy level is 0.5-8.0 μJ, and in a preferred embodiment, the pulse-on-demand mode Q-switched laser delivers an average pulse energy of 1.8 μJ, with a single pulse ranging from 1.3 μJ to 2.3 μJ. In one embodiment, the pulse width is in the range of 3ns to 1ms, preferably, in the range of 5ns to 9ns. However, those of ordinary skill in the art will know that in order to achieve the desired target accuracy and / or effect, any high-power laser with appropriate high-energy pulses and pulse frequencies that exists now or is developed later will be suitable for the present invention.

[0035] In one embodiment, in order to deliver pulsed energy from a focused energy device to a target particle or its surrounding area, a tight action region is required, which is important to minimize the potential impact of delivering energy to particles outside the target particle or region, otherwise unselected or non-target particles. For example, focused energy devices such as the Explorer XP 355-1 Q-switched laser are able to provide <4% rms, providing high pulse-to-pulse stability when fired at regular uniform intervals. However, for flow cytometry analysis and action systems where particles or cells enter the action region at non-uniform intervals, additional measures are taken to deliver uniform pulsed energy to only impact the target particles or cells or their surrounding areas. Such measures include matching laser performance parameters, such as pulse length and peak power level, so that the system of the present invention can achieve the desired target accuracy (i.e., in one embodiment, a photodamage or kill rate of 95% or more hits on the target particles). In addition, the on-demand pulse operation and performance of the laser are further adjusted to provide extremely high pulse-to-pulse stability when fired at non-uniform intervals, greatly reducing the spatial variability of the area affected by the pulses. Thus, by reducing pulse-to-pulse variability in a focused energy device, unintended motion, damage or destruction of non-target particles or cells is greatly reduced, achieving, for example, 85% or greater viability, of non-target particles or cells.

[0036] The focused energy device can be configured to damage, alter, disable, kill or destroy targeted or unwanted particles or cells in the sample fluid with pulses, or activate one of several mechanisms in the particles or cells, resulting in cell damage or death. However, depending on the desired configuration, the targeted or selected particles can be desired objects, in which case the focused energy device is not activated or triggered, or the targeted or selected particles can be unwanted particles, in which case the focused energy device is activated to act on the particles, such as to damage, alter, disable, kill or destroy the selected, unwanted particles. However, these are not the only embodiments, and various embodiments are discussed further below.

[0037] In one embodiment, when the selected particles or cells are damaged, altered, disabled, killed or destroyed by the focused energy device, the particles or cells continue to flow through the channel and flow into the container along with any non-target particles or cells. Any sheath or buffer also passes through the output channel in laminar flow.

[0038] Thus, in one embodiment, the present methods and systems are capable of producing an identified product of particles or cells in a container, including a high survival rate of non-target or desired particles or cells, and a high percentage of target particles or cells that are photodamaged, altered, disabled, destroyed or dead.

[0039] Sorting methods

[0040] On the other hand, a method for identifying subpopulations in a mixture of particles is provided. The mixture of particles may contain two or more subpopulations. These subpopulations can be distinguished from each other by one feature or multiple features. The method includes measuring a feature of a first group of particles in the mixture of particles. The method also includes distinguishing the measured feature values ​​of the subpopulations in the set of particles. On the one hand, the values ​​of the measured feature can be distinguished for two or more subpopulations in the set of particles. On the other hand, the values ​​of the measured feature can be distinguished for each subpopulation in the set of particles. The method also includes assigning a critical value range to two or more subpopulations based on the distinguishable values ​​of the measured feature. The assignment of the critical range can be performed by identifying (or classifying) subpopulations related to other subpopulations based on the distinguishable values ​​of the measured feature of each subpopulation being compared. After assigning the critical value range, the method also includes measuring a second group of particles in the mixture of particles and distinguishing the measured feature values ​​of the subpopulations in the second group of particles. Then the measured feature values ​​of the sub-particle group distinguished for the second group of particles are used to update the critical range of values ​​assigned to the sub-particle group. The method also includes using the updated critical value range to classify the particles in the mixture according to the subpopulation classification.

[0041] In one aspect, the method of the present invention is performed on a mixture of particles to identify, distinguish, differentiate or classify the particles in the mixture. The mixture of particles may contain two or more subpopulations. These subpopulations may be distinguished from each other by one feature or multiple features. These particles may be cells. These cells may be of a single type or a mixture of multiple cell types, such as a mixture of cells derived from tissues or present in blood or blood fractions. Preferably, the cells may be sperm cells. Sperm cells may be from any type of animal. More preferably, sperm cells are from non-human mammals. Further preferably, sperm cells are from livestock. Livestock include, but are not limited to, poultry, ungulates, ruminants. Further more preferably, sperm cells are sheep, goats, horses, pigs or cattle. Bovine sperm cells may include sperm cells of domestic cattle, bison, African buffalo, water buffalo or yaks.

[0042] In another aspect, the characteristics of a particle or a group of particles are detected and measured. The characteristic can be any aspect of a particle that can be detected. The characteristics can include aspects such as size, shape, permeability, integrity, content, vitality, mortality, metabolic state, and surface markers. In another aspect, the characteristic can include DNA content. DNA content, particularly in haploid cells, such as sperm cells or eggs, is an indicator of the presence of one of the two sex chromosomes. The X chromosome is larger than the Y chromosome, so haploid cells containing X chromosomes contain more total DNA than gametes containing Y chromosomes. Characteristics can be detected using various methods known to those skilled in the art, including light scattering, interference, absorbance, or fluorescence detection. For example, the characteristics of the cell DNA content can be detected by staining the cell with a DNA-binding fluorescent dye and detecting the fluorescence of the dye. In some aspects, the characteristics of the particles are measured using a flow cytometer. The use of flow cytometers for detecting and measuring characteristics is known in the art, but in general, flow cytometers operate by converting photon incident with a detector into an electrical signal. Therefore, when exposed to an electromagnetic radiation source, the electrical signal from the detector corresponds to the electromagnetic radiation from the particle (in the form of photons emitted or scattered by the particle). Typically, fluorescent molecules are used to specifically mark specific aspects of particles, and the cytometer detects photons of specific wavelengths emitted by particles in response to exposure to specific electromagnetic radiation. Electromagnetic radiation indicating characteristics from particles is detected by detectors in the cytometer device. The detector is typically a photomultiplier tube (PMT) or an avalanche photodiode (PMT), as described in U.S. Pat. No. 10,069,027. The detector converts photons into electrical signals, which can be transmitted to the processing components of the cytometer, thereby assigning values ​​to the detected features. Distinguishing the values ​​of the measured features based on the evaluation derived from the cytometer may further involve classifying and comparing the values ​​of the particles or particle groups to clarify the features of the subpopulations in the particle mixture. For example, it may be achieved by determining whether the values ​​represented graphically form a group or peak. Then, the identity of the subpopulation may be determined based on the relationship of any group or peak clarified. In some aspects, it may be distinguished by finding discrete event groupings based on one or more features measured. The measured features may be represented in a feature space. For example, a feature space may be a graphical representation of one or more measured features of a particle. In an exemplary embodiment, the feature space can be a histogram of fluorescence intensity, or a scatter plot of fluorescence intensity and total fluorescence, where the fluorescence indicates the amount of DNA present in the cell. In some aspects, the distinction can be a distinction between two subpopulations, a distinction between a single subpopulation among several subpopulations, or a distinction between a subpopulation within a continuum.

[0043] In a further aspect, event data of particles are generated by a flow cytometer. Event data is usually detected based on the electromagnetic radiation generated by the emission or scattering of particles when exposed to an electromagnetic radiation source. The electromagnetic radiation is converted into an electrical signal by a detector (usually a PMT or APD). Then, the electrical signal is converted into data by a processor component of the cytometer. In certain embodiments, the event data can be fluorescence intensity (i.e., the highest fluorescence value of a single event corresponding to a single particle) and / or total fluorescence (i.e., the integrated fluorescence value of the entire event corresponding to a single particle).

[0044] In certain aspects, the methods of the present invention involve discrete actions performed on specific groups of particles. For example, measuring a characteristic of particles in a mixture of particles and distinguishing the values ​​of the measured characteristic can be performed on a first group of particles from a mixture of particles, and subsequently, measuring a characteristic of particles in a mixture of particles and distinguishing the values ​​of the measured characteristic can be performed on a second group of particles from the same mixture of particles.

[0045] On the other hand, the method involves assigning a critical value range to each subpopulation based on previously distinguished values ​​of the measured feature. In some aspects, the assignment of critical ranges is used to determine a specific subpopulation in a particle mixture. For example, a critical range can be specified for a group of particles to identify this group of particles based on the relationship between the peak formed by the values ​​of these particles represented graphically and the peak values ​​of other particle values. The assignment of critical value ranges can include the calculation of discriminant values. The discriminant value can be a threshold discriminant value (a single value or two values) or a gate threshold (for example, three or more values). The calculation can use one or more values ​​of the parameter set. In some embodiments, the calculation of the discriminant value involves identifying at least one critical value of the first feature from the distribution of event data of the first group of particles; obtaining the value of the second feature associated with the critical value of the first feature; generating a distribution of values ​​for the second feature; identifying at least one critical value of the second feature from the distribution of the values ​​of the second feature; co-locating the critical values ​​of the first and second features in the feature space; and determining at least one discriminant value from the critical values ​​co-located in the feature space. In a further aspect, the calculation of the discriminant value can further involve adjusting the critical value to include event data from the second group of particles. In a further aspect, the critical value includes one or more peak characteristic values ​​of the first subpopulation and one or more peak characteristic values ​​of the second subpopulation. The peak characteristic values ​​include peak fluorescence intensity and peak total fluorescence. In a further aspect, the critical value can also be a valley between the peak characteristic values ​​of the first and second subpopulations.

[0046] On the other hand, the method may include measuring a characteristic of each subpopulation in the second group of particles in the particle mixture and distinguishing the value of the measured characteristic, which occurs after the critical value range is assigned to the first group of particles. The second group of particles is measured and distinguished to determine whether the critical value range assigned to the subpopulation is accurate, and if it is inaccurate, the value used for assignment is updated. In some aspects, based on the distinction of the measured characteristic values ​​of each subpopulation in the second group of particles, the value of the second group of particles is compared with the value of the first group of particles, and if they are not different, the critical range based on the value assigned to the first group of particles is used. If there is a difference between the values ​​of the first group of particles and the second group of particles, the critical value range is updated using the value of the second group of particles. Typically, a certain number of particles, or the number of particles that can be measured in a certain time, determines the composition of each group of particles, thereby determining the frequency of updating the critical value range. The number of particles composed of the set and / or time value is part of the parameter set.

[0047] In another aspect, classifications are assigned to individual particles in a mixture of particles based on a specified critical range. This can be based on a critical range that is updated or not updated due to a new set of particles being measured. This classification of particles is derived from the relationship between subpopulations. Based on this critical value assignment, individual particles in the mixture can be classified as belonging to a subpopulation, while these classified particles are part of a second set of particles used to update the critical range.

[0048] On the other hand, the method of the present invention provides a method for identifying a subpopulation in a particle mixture by defining a gate for distinguishing one or more particle subpopulations in the mixture. When the particles in the mixture are processed by a flow cytometer, the particles are detected by the detection device of the flow cytometer, thereby generating a detection event. Each detector event can be a collection of values, including, for example, values ​​of the intensity and area of ​​the signal. The event data is processed by a computer. The event data is accumulated until a sufficient number of event data are compiled to allow one or more discriminant values ​​to be calculated. Before calculating any discriminant value, a parameter set is provided, i.e., a set of values ​​that affect certain details of all discriminant value calculations. When enough event data is compiled, the discriminant value is calculated and the process is repeated with a new set of event data. Each set of event data is a previous set of event data to determine whether the discriminant value calculated by a subsequent set of event data is significantly different from the discriminant value calculated from the previous set of event data. If the difference is not significant, the available data obtained by the previous identification value calculation can be averaged to improve the consistency of position and shape between the continuously calculated identification values. If the discriminant values ​​are significantly different, the data of the previously calculated discriminant value are discarded.

[0049] In certain aspects of the method, a mixture of particles is provided. The composition of the mixture may vary depending on the application. The mixture may include particles and a medium. In certain embodiments, the mixture of particles is a sperm sample. The sperm cells may contain raw ejaculate components, other medium and / or thickeners. The particle mixture is delivered to a flow cytometer system.

[0050] In certain aspects, a parameter set is provided. The parameter set is provided to a processor for calculating one or more discriminant values. The parameter set includes a set of values ​​that affect certain details of all discriminant value calculations. In certain aspects, the parameter set includes values ​​for calculating thresholds. In other aspects, the parameter set includes values ​​for calculating one or more points or vertices of a gate. In further aspects, the parameter set includes a frequency at which discriminant value calculations occur and / or an amount of event data that triggers the discriminant value calculations.

[0051] In some aspects, the characteristics of particles in the mixture are detected. The detection is performed by a cytometer. The cytometer detects the characteristics of a single particle by detecting various aspects of the particle. The detection of a single cell characteristic is considered to be a single event. In some cases, the particles may arrive at the detector in the form of a particle beam, so for aggregated cells, the detection will be performed as a single event. These aspects may include the size and shape of the particles (e.g., by detecting light scattering) and the composition of the particles (e.g., by detecting the presence of fluorescent markers on the surface of the particles or within the particles). In some aspects, the detected particle aspects may be fluorescent dyes. In further aspects, the particles may be cells, and the dye may be a fluorescent DNA dye. One or more characteristics may be detected for each particle, and a single event may include one or more detected characteristics. The detector of the cytometer generates a signal corresponding to the detected characteristics of each event and provides the signal to a processor (e.g., a CPU).

[0052] At step 306, the compiled event data is distributed, and the occurrence of data distribution may depend on one or more values ​​of the parameter set - the frequency or number of events. In some aspects, the distribution can be to sort the events by event area. In other aspects, the distribution can be to sort the events by intensity. In further aspects, the distribution can be to sort the events by a classification of intensity and area.

[0053] In another aspect, the classification of events is provided to a sorting device that implements a sort / no-sort decision.

[0054] Compare the discriminant value calculated for a compiled set of events with the discriminant value calculated for the previous set of events. [Repeat steps 303, 307 to generate new discriminant values ​​based on the newly compiled set of event data]. [If different, discard the previous; if no difference, average the event data]. [This can be done by directly comparing the calculated discriminant values, or by comparing the data in the two compiled data sets].

[0055] In another aspect, the method of the present invention can be carried out as follows:

[0056] 1. Get the parameter set.

[0057] 2. Get the frequency of gate calculation from the parameter set.

[0058] 3. Get the maximum number of events from the parameter set whose data will be used directly in the gate calculation (Note: Each event has an "intensity" which is a measure of the maximum instantaneous fluorescence, and an "area" which is the integral of the intensity over the duration of the event.)

[0059] 4. As the detector generates events, each event is added to the event data set used to calculate the gate. If the maximum number of events already exists, the new event is added and the oldest event is ignored in subsequent gate determinations.

[0060] 5. According to step 2, periodically determine whether enough events have accumulated to attempt to calculate a gate, and if so, continue to step 6.

[0061] 6. Sort the accumulated events by region. Save the median value of the event region in this dataset.

[0062] 7. Sort the accumulated events by intensity and then by area. Save the median intensity of the events in this dataset.

[0063] 8. Determine if the characteristics of the event data have changed significantly since the last calculation by comparing the median values ​​of intensity and area to the corresponding median values ​​in the last calculated data set. (If no previous median values ​​were available, proceed as if there were no significant changes.)

[0064] 9. If the characteristics of the event data have changed significantly (from step 8), first save and then clear the data associated with the previous gates so that it does not affect the gates generated by the current dataset. For example, "forget" old data that may have been averaged with the current data in steps 32, 57, and all other similar steps.

[0065] 10. Create a histogram of the intensity values ​​of all events. The values ​​placed in the histogram bins range from 0 to 65535, and they are sorted into 512 equally sized bins that together cover the entire range.

[0066] 11. Calculate 25% from the left edge of the histogram and 15% from the right edge of the histogram. For 512 bins, these numbers correspond to 128 and 435.

[0067] 12. Find the three highest peaks in the histogram region bounded by the limits determined in step 11. A peak in the histogram is one where the number of counts is greater than the counts in the bars on either side.

[0068] 13. Sort the bin numbers of the histogram peaks in step 12 in descending order by bin height (number of events in each bin).

[0069] 14. Check all peaks within the histogram region outside the limits determined in step 11. If any of them are higher than the second highest peak determined in step 13, flag a gate calculation error and wait for the setup to become available for new event data.

[0070] 15. There will be some expected species-dependent range of intensity differences between the male and female peaks. Call this the desired intensity difference. Determine the intensity differences between the three highest peaks from step 13 and take pairs of two peaks. There are three pairs: peaks 1 and 2, peaks 1 and 3, and peaks 2 and 3.

[0071] 16. If no pair of peaks show the required intensity difference, mark the calculation error and wait for a newly available event data set.

[0072] 17. If the difference between only one pair in step 15 is the required intensity difference, then the pair of histogram peaks represents the male and female peaks. If the third peak that is not in the pair is the highest of the three peaks in step 12, then mark the calculation error and wait for a new available event data set.

[0073] 18. If step 17 did not identify male and female histogram peaks, determine if both pairs of the highest peaks show the desired intensity difference. If two pairs do, then all three peaks are equally spaced by bin number (case A), or the two highest peaks are very close together (one is a sharp peak, both represent the same male or female peak), and the third represents a peak of the opposite sex (case B). Check the spacing to determine if case A or B is correct.

[0074] 19. If you encounter case A in step 18, the central peak is either a male peak or a female peak. Of the left and right peaks, the height of one must be less than two-thirds of the height of the other peak. If this is not the case, mark a calculation error and wait for a new available event data set. If this is correct, the higher of the left and right peaks is the opposite sex peak; the higher of the central peak and the two outer peaks are the male and female peaks. The peak on the left may be a male peak or a female peak, depending on the kind of sample source.

[0075] 20. If situation B in step 18 is encountered, the histogram peak among the two pairs of peaks that exhibit the required intensity difference is selected as one of the related (male or female) histogram peaks, and the other two histogram peaks are regarded as another related peak and further processed to determine which is the peak of the peak and which is the main peak (see step 21).

[0076] 21. If you encounter case B in step 18, find the width of the two histogram peaks representing the remaining peak of interest and the peak of this peak. Calculate the width by finding the bar with the lowest point (trough) between them, then move outward from the trough to find the first bar in both directions that is lower in height than the trough. The number of bins from the trough bin to the first lower bin is the width of the peak in that direction. Take the narrower of the two peaks as the sharp peak and the wider one as the second peak of interest (male or female peak).

[0077] 22. If the two peaks identified as male and female are not the two highest peaks, mark the calculation error and wait for a new available event data set. If so, designate the rightmost of the two peaks as the female intensity peak and the leftmost peak as the male intensity peak, or vice versa, depending on the species from which the sample originated.

[0078] 23. Identify the lowest histogram bar between the identified male and female peaks as the valley bottom.

[0079] 24. Calculate the peak-to-valley ratio (P2T ratio) as ((P1+P2) / 2) / T, where P1 and P2 are the counts in the peak bins of the male and female intensity histograms, respectively, and T is the counts in the trough bin. In addition, retrieve the peak-to-valley threshold (P2T threshold) from the parameter set. In addition, calculate the peak dominance value as (P2-T) / (P1-T), and retrieve the minimum peak dominance value from the parameter set.

[0080] 25. Compare the P2T ratio calculated in step 24 with the P2T threshold. If the peak-to-valley ratio is less than the threshold, mark an error and wait for a new available data set. At the same time, compare the peak dominance value and the minimum peak dominance value. If the peak dominance value is less than the minimum peak dominance value, mark an error and wait for a new available data set.

[0081] 26. For the male and female peaks, count the number of bins from the valley to the peak of the histogram.

[0082] 27. For the male and female histogram peaks, calculate the width of the peak in the bin as twice the corresponding value found in step 26, plus 1.

[0083] 28. For the male and female histogram peaks, calculate the data for the histogram bins using the curve fit point (x, y = binnumber, count). Calculate this as the peak width in step 27, unless the peak width is greater than 9, otherwise set the number of stored bins to 9.

[0084] 29. For the male and female histogram peaks, find a more precise male and female peak location by applying a polynomial best fit algorithm to the bins around the peaks, using the number of bins calculated in step 28 as the number of points for the curve fit; then apply the Newton-Raphson algorithm to find the intensity corresponding to the zero slope position in the fitted polynomial. This more precise peak estimate is a non-integer number of bins with a decimal point.

[0085] 30. Starting from step 23, apply the polynomial best fit algorithm to the five bins centered at the trough bin. Then apply the Newton-Raphson algorithm to find the intensity corresponding to the zero slope position in the fitted polynomial. This more accurate valley bottom estimate is a non-integer bin number with a decimal point.

[0086] 31. If the difference between the bin estimate in step 30 and the bin quantity initially found in step 23 is greater than 1.0, adjust the refined estimate to the bin value in step 23 plus or minus 1.0, depending on whether the refined estimate in step 30 is greater than or less than the bin value in step 23.

[0087] 32. Calculate the average of the bin values ​​for step 31 and the bin values ​​for step 31 obtained in the last five event data sets. Use this average as the trough value obtained from this data set.

[0088] 33. Multiply the width of the intensity value range (65536) by 0.003. The product is 196.608.

[0089] 34. Convert the male and female histogram peak positions (as bin numbers) to intensity values ​​by multiplying the exact bin number by the bin width. For example, if the female histogram peak position in step 29 is 257.358, and there are 512 histogram bins in the range 0 to 65535, then the bin width is 128, and the intensity peak in the range 0 to 65535 is obtained by multiplying 257.358×128=32942.6.

[0090] 35. Using the result of the calculation in step 33 and the product of the peak positions of the male and female histograms converted to intensity values ​​from step 34, find a subset of events whose intensity values ​​are within the peak distance given by the product in step 33. Do this separately for the male peak and the female peak. For example, if the female intensity peak is 32942.6, find a subset of events in the event data set whose intensities are in the range [32942.6–196.608, 32942.6+196.608], that is, events with intensities between 32746 and 33139. Do the same for a set of events whose intensities are close to the male histogram peak.

[0091] 36. Perform the following steps 37 to 54 for each male and female peak to find the "area" coordinates of the peaks to be gated in the (intensity, area) plane. These steps use only the event data from the events selected in step 35.

[0092] 37. Create an array of area values ​​in the event selected in step 35. Sort in ascending order.

[0093] 38. Calculate a number equal to 75% of the number of area values ​​in the array from step 37. In this process description, this will be referred to as numValuesIn75Pct.

[0094] 39. Calculate a number that is equal to 99% of the number of area values ​​in the array from step 37. In this process description, this will be called largestRetainedIndex.

[0095] 40. In subsequent steps, ignore all area values ​​greater than the 99th percentile value, that is, all values ​​in the array whose array index is greater than largestRetainedIndex (the index calculated in step 39).

[0096] 41. Find the area value that marks the 75th percentile of area values, using numValuesIn75pct from step 38. In this process description, this value will be referred to as "lowest75pctMax".

[0097] 42. Find the minimum area value in the array from step 37. In this process description, this value will be referred to as "lowest75pctMin".

[0098] 43. Calculate lowest75pctMax–lowest75pctMin, which represents the range of area values ​​that are not ignored, as described in step 40. In this process description, this difference is referred to as “areaRange”.

[0099] 44. Calculate the variable width of the "area" histogram bars as: areaBinWidth = areaRange / 512.

[0100] 45. Calculate the value (largestRetainedArea / areaBinWidth)+1, which will be referred to as "numberOfAreaBins" in this process description.

[0101] 46. ​​Look up the value at array index "largestRetainedIndex" in step 37. In this process description, this value will be referred to as "largestRetainedArea".

[0102] 47. If there are few selected events to view, then:

[0103] largestRetainedIndex<(6x numberOfAreaBins from step 45)

[0104] then recalculate areaBinWidth according to:

[0105] numberOfAreaBins=(largestRetainedIndex / 6)+1

[0106] If(numberOfAreaBins<4)

[0107] numberOfAreaBins = 4;

[0108] areaBinWidth=largestRetainedArea / (numberOfAreaBins–1)

[0109] 48. Create a histogram of the area values ​​in the array of step 37, ignoring all values ​​in the array whose array index is greater than the maximum retained index mentioned in step 40.

[0110] 49. Create another histogram of area values ​​using the same number of bins, smoothed by averaging the value counts in each bin from step 48 with the value counts in the seven bins on either side (from step 48), i.e., adding the 15 values ​​together for the average.

[0111] 50. Find the bin with the most counts in the smoothed area histogram in step 49.

[0112] 51. Multiply the counts in the bar from step 50 by 0.2.

[0113] 52. Moving leftward (smaller array indices) from the bin with the most counts, find the first bin with a count less than the product calculated in step 51.

[0114] 53. Moving rightward (larger array index) from the bin with the most counts, find the first bin with a count less than the product calculated in step 51.

[0115] 54. Use a least squares polynomial curve fitting algorithm to find a 10th-degree polynomial that approximates the smoothed histogram area between the bins found in steps 52 and 53. Using the Newton-Raphson algorithm, find the position where the slope of the polynomial is zero and use this position as the area value of the peak.

[0116] 55. Once the (intensity, area) positions of the male and female peaks are known, calculate the position of the point in the (intensity, area) coordinate space that is on a straight line through the male and female peaks and has the intensity coordinates of the trough calculated in step 32. This point is the trough bin position in the (intensity, area) coordinate space.

[0117] 56. If the area coordinate of the male peak is higher than the area coordinate of the female peak, or lower than the area coordinate of the female peak, then mark an error and wait for the next event data set, depending on the type of sample source.

[0118] 57. For each female peak intensity coordinate in the most recent thirty event data sets, calculate the female peak intensity coordinate value of step 29 that has been adjusted for drift in the electronics that has occurred since the event data sets were processed. If fewer than thirty previous values ​​are available, simply average the available values. In this process description, the average value will be referred to as the "average female peak intensity." (The process of calculating drift-adjusted coordinates is described as a set of steps that are enumerated with letters rather than numbers.)

[0119] 58. Calculate the value of the male intensity coordinate in step 29 that has been adjusted for drift in the electronics that has occurred since the event data set was processed, for each male peak intensity coordinate in the thirty most recent event data sets. Average the drift-adjusted male intensity peak coordinate with the male intensity coordinates in step 29 for processing this event data set. If fewer than thirty previous values ​​are available, then average only the available values. In this process description, the average will be referred to as the "average male peak intensity." (The process of calculating drift-adjusted coordinates is described as a set of steps that are enumerated with letters rather than numbers.)

[0120] 59. Compute the bin intensity coordinate value of step 32, adjusted for drift in the electronics that has occurred since the event data set was processed, for each trough intensity coordinate of the last thirty event data sets. Average the drift-adjusted bin intensity coordinates with the bin intensity coordinates from step 32 for processing this data set. If fewer than thirty previous values ​​are available, average only the available values. In this process description, the average is referred to as the "average trough intensity." (The process of calculating drift-adjusted coordinates is described as a set of steps that are enumerated with letters rather than numbers.)

[0121] 60. Calculate the difference between the bin intensity in step 32 and the peak intensity of the gate drawn in step 29. Call this distance the peak-to-valley intensity distance.

[0122] 61. Average the peak-to-valley intensity distance of step 60 with the peak-to-valley intensity distances of the most recent 30 event data sets found in step 60. If fewer than 30 previous values ​​are available, then only the available values ​​are averaged. In this process description, the average is referred to as the "average peak-to-valley intensity distance".

[0123] 62. Create a polygonal gate with eight vertices in the (intensity, area) coordinate plane based on the locations of the vertices of interest. Start this stage of the process by getting the encroachmentPct value from the parameter set. Divide by 100.0 and then subtract the result from 1.0 to get the "encroachment multiplier" used in subsequent steps.

[0124] 63. Multiply the average peak-to-valley intensity distance by the erosion multiplier from step 62. In this process description, the product is referred to as the "Encroachment Adjustment."

[0125] 64. A value called "Erosion Strength" is calculated for use in this process description. The two vertices lie on a vertical line drawn through the erosion strength. If the gate was created relative to the male peak, the erosion strength is the average male peak intensity plus the erosion adjustment. For the male peak, the erosion strength is the intensity coordinate of the gate vertex or vertices that are farthest to the right on the intensity axis in the (intensity, area) plane. If the gate was created relative to the female peak, the erosion strength is the average female peak intensity minus the erosion adjustment. For the female peak, the erosion strength is the intensity coordinate of the gate vertex or vertices that are farthest to the left on the intensity axis in the (intensity, area) plane.

[0126] 65. Retrieve the peakFalloffPct value from the parameter set. This is a decimal value in the range 0.0 to 100.0. Currently 65.0 is used, but eventually some other value may be determined to be optimal. Divide peakFalloffPct by 100.0 and then multiply it by the height of the highest bar in the intensity histogram for the associated peak (step 10). The associated peak is either the male intensity peak or the female intensity peak, depending on which peak the gate is calculated relative to. In this process description, the product of the multiplication will be referred to as "falloffValue".

[0127] 66. Starting from the highest intensity histogram bar in the correlation peak, move away from the bin in step 23 until the height of the currently examined bar is less than or equal to falloffValue. Calculate the intensity value corresponding to the center of the bar as [bar index + 0.5] × [bar width]. In this process description, this value is called the "falloff intensity".

[0128] 67. Average the decay intensity from step 66 with the decay intensities of the last thirty event data sets. If fewer than 30 values ​​are available previously, only the available values ​​are averaged. In this process description, the average is referred to as the "average decay intensity".

[0129] 68. Take the average decay intensity as the initial value of the intensity coordinate in the (intensity, area) plane for the two gate vertices that are on the rightmost side in a gate drawn relative to the female peak or on the leftmost side in a gate drawn relative to the male peak.

[0130] 69. Retrieve the "octExtension" value from the parameter set. octExtension is a decimal value greater than 0.0. Currently a value of 0.6 is used, but this may change in the future. Multiply octExtension by the average peak-to-valley intensity distance from step 61. In this process description, this value will be referred to as the "extension distance."

[0131] 70. If the gate is drawn relative to the male peak, subtract the extension distance in step 69 from the initial value calculated in step 68. If drawn relative to the female peak, add the extension distance to the initial value. Use the result as the final value for the intensity coordinate in the (intensity, area) plane for the two gate vertices that are the rightmost in a gate drawn relative to the female peak or the leftmost in a gate drawn relative to the male peak.

[0132] 71. Calculate the mean invasion intensity and mean female peak intensity for gates plotted relative to the female peak, or the mean male peak intensity for gates plotted relative to the male peak.

[0133] 72. Retrieve the "horizadjmultipler1" value from the parameter set. Multiply it by the mean peak-to-valley intensity distance from step 61. Add the product to the value calculated in step 71. For gates relative to the female peak, the result is used as the intensity coordinates of the two vertices on a vertical line in the (intensity, area) plane that lies to the right of the encroaching intensity. For gates relative to the male peak, the line lies to the left of the encroaching intensity. A value of 0.0 is currently used for horizadjmultipler1, but this may change.

[0134] 73. Calculate the average of the mean decay intensity and the mean female peak intensity for gates drawn relative to the female peak, or the mean male peak intensity for gates drawn relative to the male peak.

[0135] 74. Retrieve the "horizadjmultipler2" value from the parameter set. Multiply it by the mean peak-to-valley intensity distance from step 61. Add the product to the value calculated in step 73. For a gate relative to the female peak, the result is used as the intensity coordinates in the (intensity, area) plane of the two vertices on a vertical line to the right of the line calculated in step 72 and to the left of a line drawn through the decaying intensities. For a gate relative to the male peak, this line is to the left of the line calculated in step 72 and to the right of a line drawn through the decaying intensities. A value of 0.0 is currently used for horizadjmultipler2, but this may change.

[0136] 75. Calculate the difference between the trough "area" coordinate from step 55 and the peak area coordinate of the peak for which you drew the gate in step 29. Call this distance the peak-to-valley area distance.

[0137] 76. Average the peak-to-valley area distance in step 75 with the peak-to-valley area distances of the most recent 30 event data sets found in step 75. If fewer than 30 previous values ​​are available, then only the available values ​​are averaged. In this process description, the average is referred to as the "average peak-to-valley area distance".

[0138] 77. Consider an "axis" that passes through the (intensity, area) male and female peaks. Two gate vertices lie on the line whose intensity was calculated in step 72, one above the axis and one below the axis. Calculate the intersection point (intensity, area) where the axis and the vertical line from step 72 intersect. Then calculate the "area" coordinates of the vertices on the vertical line and above the axis by retrieving the value of "AreaAboverMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 76, and adding the result to the area coordinate of the intersection point.

[0139] 78. Calculate the "area" coordinates of the vertices below the axis in step 72 by retrieving the value of "areaBelowMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 76, then subtracting the result from the area coordinates of the intersection points.

[0140] 79. Two vertices lie on a vertical line drawn through the encroachment point, one above the axis line and one below the axis line. Compute the point where the vertical line intersects the axis line at the encroachment strength; in this process description, this point will be referred to as the "encroachment point." Compute the "area" coordinate of the vertex with the erosion strength and above the axis line by retrieving the value of AboutHeightPctLeft from the parameter set, dividing it by 100.0, and then multiplying it by the length of the vertical line segment from the vertex to the axis line calculated in step 77; then add this value to the "area" coordinate of the encroachment point.

[0141] 80. Starting at the point where the vertical line at the erosion strength intersects the axis line (calculated in step 79), calculate the vertex of step 79 that lies on the vertical line at the erosion strength but below the axis line by retrieving the value of belowHeightPctLeft from the parameter set, dividing by 100.0, and multiplying by the length of the vertical line segment from the vertex to the axis line calculated from step 78; then subtract that value from the area value of the encroached point. This is the "area" coordinate of the point.

[0142] 81. The two vertices lie on the perpendicular line with the intensities calculated in step 74. Calculate the intersection (intensity, area) where the axis from step 74 and the perpendicular line intersect. Then calculate the "area" coordinate of the vertex on the axis by retrieving the value of "AreaAboverMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 76, then adding the result to the area coordinate of the intersection point.

[0143] 82. Calculate the "area" coordinates of the vertices on the line that are below the axis in step 74 by retrieving the value of "areaBelowMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 76, then subtracting the result from the area coordinates of the intersection points calculated in step 81.

[0144] 83. Two vertices lie on the vertical line drawn through the falloff strength, one above the axis line and one below the axis line. Compute the point where the vertical line intersects the axis line at the falloff strength; in this flow description, this point is called the "falloff point". Compute the "area" coordinate of the vertex with the falloff strength that is above the axis line by retrieving the value of AboutHeightPCTright from the parameter set, dividing it by 100.0, and multiplying it by the length of the vertical line segment from the vertex to the axis line calculated in step 81; then add this value to the "area" coordinate of the falloff point.

[0145] 84. Starting at the falloff point from step 83, calculate the vertex from step 83 that lies on the vertical line at the falloff strength, but below the axis line, by retrieving the value of belowheightptright from the parameter set, dividing by 100.0, and multiplying by the length of the vertical line segment from the vertex to the axis line calculated from step 82; then subtract that value from the area value at the falloff point, i.e. the "area" coordinate of the point below the axis.

[0146] 85. Eight gate vertices have now been calculated from the event dataset, unless an error was flagged midway through the process. Retrieve the value for the perturbation skip from the parameter set. If the number of event datasets processed since the last significant change in the data was noted in step 9 is less than or equal to the perturbation skip (including the current dataset), do not create a gate. Instead, use "No Gate" and apply the same operation to all cells flowing through the flow cytometer system until a new gate is created by processing subsequent datasets. However, if the number of event datasets processed since the last significant change is greater than the perturbation skip, apply the gate calculated from the event dataset by this procedure.

[0147] 86. For each value that needs to be averaged with similar values ​​in subsequent iterations of the process, save the value that results from processing this data set.

[0148] In another aspect, the process also includes an inventive technique for finding drift correction strength [an expansion of steps 57-59 above]

[0149] A. Create an event dataset counter variable and initialize it to zero.

[0150] B. Increment the counter each time a new event data set is processed. Increment the counter to determine if an error was flagged when processing the data and if a character change in the event data was detected in step 8.

[0151] C. When saving a value to be averaged with similar values ​​from a subsequent data set for future use in step 86, the value of the data set counter associated with the data set that produced the value is saved.

[0152] D. Calculate the average drift rate of variable X over the period of time spanned by processing multiple data sets by subtracting the earliest available saved value of X from the latest value of X, obtaining DELTA_X, and subtracting the value of the data set counter associated with the earliest value from the value of the counter associated with the most recently stored value, obtaining DELTA_T; then calculate DELTA_X / DELTA_T as the average drift.

[0153] E. To calculate the total drift adjustment for a particular recorded instance of a variable, multiply the average drift rate for the variable obtained from step D by the number of time periods that have passed since that value was recorded (i.e., the number of event data sets that have been presented). The product of the multiplications represents the drift adjustment. Add the drift-adjusted value to the previously recorded value to obtain the drift-corrected value.

[0154] In certain aspects, the present invention provides sorted particle products produced by the systems and methods described herein. These sorted products may have higher and / or more consistent quality. In certain aspects, the sorted particle products of the present invention are sex-identified sperm cell samples. The sex-identified sperm sample may be mainly sperm with X chromosomes, or it may be mainly sperm with Y chromosomes. The sex-identified semen sample may be at least 60%, at least 65%, at least 70%, at least 75%, at least 80%, at least 85%, at least 90%, at least 91%, at least 92%, at least 93%, at least 94%, at least 95%, at least 96%, at least 97%, at least 98%, or at least 99% (or any percentage in the range of 60-100%) sperm cells with X chromosomes or Y chromosomes. The sex-identified sperm cell sample may contain a salt-based culture medium formula, optionally including any number of beneficial ingredients, including but not limited to cryopreservatives (e.g., glycerol), antioxidants, fillers, proteins, vitamins, and sugars. Sexed sperm samples may be suitable for fertilization, including artificial insemination (AI) or in vitro fertilization (IVF). Example

[0155] Identification of sperm cells carrying the Y chromosome using a program-optimized gating system

[0156] For some methods of semen sexing, the mixture of cells present in the original ejaculate (along with any dead cells and cell clusters) is enriched for sperm cells that carry a live X chromosome or sperm cells that carry a Y chromosome. For dairy applications, sperm cells that carry an X chromosome are more likely to result in pregnancy and will produce female offspring. Enrichment is usually performed by droplet sorting or laser ablation. In laser ablation applications, the "kill" laser attempts to kill all cells that do not have the desired characteristics (i.e., any cells that are not sperm cells that carry a live X chromosome). In the current system, cells are labeled with a DNA-intercalating fluorescent dye, and the fluorescence emitted by each cell is detected as a measure of the total DNA content. The fluorescence signal generated by each cell (or cell cluster) is detected by an avalanche photodiode and recorded as an event. Bovine sperm cells that carry an X chromosome have approximately 3.8% more DNA than sperm cells that carry a Y chromosome, which means that sperm cells that carry an X chromosome will fluoresce more. The fluorescence intensity (peak) corresponding to sperm cells that carry an X chromosome is greater than that of sperm cells that carry a Y chromosome. The area of ​​the fluorescence signal for each event is also determined. The fluorescence signal area is the integral (i.e., sum) of all fluorescence from that event over the amount of time that the event was detected. Each event can be graphically represented on an intensity (peak) and area graph, and when graphs for many events are compiled, subpopulations of cells can be identified. For example, a subpopulation of sperm cells bearing an X chromosome will have a higher density relative to a subpopulation of sperm cells bearing a Y chromosome and as a single cell based on area (a dichromosome cell will have a higher density and / or a larger area). The kill laser is triggered to fire any event when the event data produces results that fall outside a certain area on the graph (i.e., the area containing sperm cells bearing an X chromosome). The boundaries of this area on the graph are the gates. In previous iterations of the system, the software operator created gates by clicking on locations in the peak intensity graph. The polygon describing the gate was created by selecting vertices in sequence until the polygon was complete. If the polygon was successfully described (i.e., it contained at least three vertices), the coordinates of the gate vertices were sent to the GPP and from the GPP to the DSP to control the kill laser so that the kill laser is triggered by an even number of events detected that fall outside the gate.

[0157] It has been observed that different operators draw gates in different ways, resulting in unknown amounts of variation in a particular desired outcome, whether it be increased purity or increased yield of sexed semen samples. Furthermore, previous work has demonstrated that gates attributed to a subpopulation at the beginning of a run may not attribute to the same population at the end of the run, which can be due to several factors. Thus, there are three issues that need to be addressed in previous systems: lack of control over the shape and position of gate generation, resulting in inconsistent products and suboptimal collection conditions; inadequate centering of gates around maternal peaks due to drift in event attributes during the collection process; and operator involvement to establish and maintain matching of gates to the population of sperm cells carrying X chromosomes.

[0158] To address these issues, an algorithm for drawing gates by setting gate vertices is proposed based on the known relationship between the X chromosome sperm population and the Y chromosome sperm population. Specifically, when compiling events, the gate vertices are determined based on the distribution of events by using an algorithm that performs the following operations. In terms of the performance of the algorithm, the system uses a parameter set including: the frequency of gate calculation and the maximum number of events that will be directly used for gate calculation; the peak to trough (P2T) threshold, and the configuration settings used to calculate the gate vertices. The algorithm is executed as follows:

[0159] 1. Accumulate events – Events are added one at a time from an external source (the flow cytometer) to the register by calling a static function, which in turn calls a public member function. The frequency of gate calculation and the maximum number of events that can be used directly for gate calculation are obtained from the parameter collection.

[0160] a. The events are stored in an array. If the number of events in the array is less than the default value set when the function is called, an error is returned.

[0161] b. The array of accumulated events is a circular buffer. The number of events currently in the buffer (starting at index 0) is initialized to zero and incremented each time an event is registered. Once the number of events reaches the size of the buffer, it is no longer incremented and new events added will overwrite the oldest events. The second variable provided is the index of the next event to be written. It is also initialized to zero.

[0162] c. Determine if enough events have accumulated to attempt to calculate the gate. Once enough events have accumulated, proceed with calculating the gate.

[0163] d. The accumulated events are copied into a separate array and then sorted in ascending order by peak and area values. Sorting is achieved by creating an array of pointers to the events. The median value of the area of ​​the events in this data set is saved.

[0164] e. It is also possible to add multiple events at the same time by providing an array of events. However, this is not a static function and the caller needs to call it. The caller decides whether to add them individually or in a buffer of multiple events.

[0165] f. Determine whether the characteristics of the event data have changed significantly since the last calculation by comparing the median values ​​of intensity and area with the corresponding median values ​​in the last calculated data set. If the characteristics of the event data have changed significantly, first save and then clear the data associated with the previous gate so that it does not affect the gate generated by the current data set.

[0166] 2. Calculate the gate by calling calculategate function. The following occurs in this function:

[0167] a. Create a histogram of event peaks. The values ​​placed in the histogram bins range from 0 to 65535, and they are binned into 512 equally sized bins that together cover the entire range. There are 512 bins, and each bin is 128 wide.

[0168] b. Find the highest peak in the histogram, if any, and determine which meets the criteria for male and female peaks. Find the three highest peaks and sort them by height. Three were determined instead of two to account for the situation where the two highest apparent peaks are not actually one peak, but rather are male and female peaks or peaks flanking the male peak.

[0169] 1. Calculate the indices that allow entry into the leftmost and rightmost peak position histograms. These indices are the 25% left edge of the histogram and the 15% right edge of the histogram.

[0170] 2. Determine the three highest peaks within the allowed region based on the bar numbers. A "peak" is a bar graph bar with lower values ​​on both the left and right bars. If two or more adjacent bars have the same value, and both bars are lower, they cannot be identified as peaks at this time and should be considered an error.

[0171] 3. Sort the bin numbers of the histogram peaks from the previous step in descending order by bin height (number of events in each bin).

[0172] 4. Check the peaks in the portion of the histogram that are outside the allowed range of the leftmost and rightmost allowed peak positions. If any peak found outside the allowed range is higher than the second highest peak of the three highest peaks found within the allowed range, an error is returned because the male and female peaks are required to be the two highest peaks and must be within the cutoff limits on the left and right.

[0173] 5. If the male and female peaks are correctly identified, the distance between them is always between 2.8% and 4.8%. The algorithm uses this information to process the three highest peaks, if any, and determine which one is the male peak and the female peak. If only two peaks are found within the allowed range, the distance between them must be correct, otherwise the calculation error is marked and a new available event data set is awaited.

[0174] 6. If three peaks are found within the allowed range, determine which of the two near peaks is the spike of the other by finding the lowest point between the two near peaks, then moving away from the lowest point in both directions to find the other side of the two near peaks at the same height as the bin. The wide side is the "real" peak and the narrow side is the "spike".

[0175] 7. If exactly two of the three highest peaks are 2.8% to 4.8% apart, then these two peaks are candidates for the male and female peaks. Check the third peak to see if it is one of the other two peaks because it is very close to one of them. If it is a spike, ignore it. If it is not a spike and is higher than any of the other two peaks, mark the calculation as wrong and wait for a new available event data set.

[0176] 8. If the spacing between these three peaks is equal and about 3.8%, then

[0177] i. If the leftmost peak is less than 2 / 3 of the height of the rightmost peak, then the middle and rightmost peaks are candidates for male and female peaks; or

[0178] ii. If the rightmost peak is less than 2 / 3 of the height of the leftmost peak, then the middlemost and leftmost peaks are candidates for the male and female peaks.

[0179] iii. If neither (i) nor (ii) is satisfied, an error is returned.

[0180] 9. If the distances between the three highest peaks are not all within about 3.8%, flag a calculation error and wait for a new available event dataset.

[0181] 10. After the candidate peaks for male and female peaks are selected by any of the above methods, they are double checked to confirm that they are the two highest peaks, or the other, the highest peak is the peak of one of the other two peaks. If this is not the case, an error is generated. Otherwise, the candidate peaks are considered to be male and female peaks.

[0182] c. Determine the x-coordinate (intensity) value of the peak. Calculate the exact x-coordinate peak for the identified male and female histogram peaks to be more accurate than their bar identifiers.

[0183] 1. Apply a least squares polygon fit algorithm to determine the exact histogram cell for the center of each peak. The peaks in the histogram are very sharp. Only a few bins are used to calculate the center of the peaks. The leftmost and rightmost bins to use are determined by calculating the distance from peak to valley between the two peaks (the number of histogram bins), and using the distance from the peak center in both directions to determine the limits of the histogram peak bins that are input to the curve fit algorithm.

[0184] 2. The degree of the polynomial used for fitting is set to the number of bars entered into the algorithm. After fitting the polynomial, apply the Newton-Raphson algorithm to find the value of "x" near the estimated peak where the slope of the polynomial is zero, with an accuracy of at least 1 part in 65536. The result of the Newton-Raphson algorithm is used as the center of the peak.

[0185] 3. If a NaN (not a number) occurs in a floating point processor while fitting a polynomial curve, an error is returned but execution continues. If the Newton-Raphson algorithm cannot determine the center of the peak within 1000 iterations (usually only about 10 iterations are required), an error is returned but execution continues.

[0186] 4. The x-coordinate of each peak is calculated as the exact histogram bin of each peak center provided by the polynomial fitting algorithm multiplied by the histogram bin width.

[0187] 5. The female peak is considered to be the rightmost of the two identified peaks, ie the peak with the larger event peak "x" coordinate.

[0188] d. Identify the bins that contain the trough. The peak-to-trough ratio (P2T) is calculated and compared to the P2T threshold, otherwise set to the default value of 1.2. If the P2T ratio in the histogram is too low, an error is generated.

[0189] 1. Calculate the P2T ratio as ((P1+P2) / 2) / T, where P1 and P2 are the counts in the peak bins of the male and female intensity histograms, respectively, and T is the count in the trough bin.

[0190] 2. Calculate the peak dominance value as (P2-T) / (P1-T) and retrieve the peak dominance threshold from the parameter set.

[0191] 3. Compare the calculated P2T ratio with the P2T threshold. If the peak-to-valley ratio is less than the threshold, mark an error and wait for a new available dataset.

[0192] 4. Compare the calculated peak dominance value with the dominant peak dominance threshold. If the calculated peak dominance value is less than the threshold, mark an error and wait for a new available dataset.

[0193] e. Determine the bin height as the histogram height of the lowest point between the two peaks in the histogram.

[0194] f. Create a histogram of the event area values ​​for both male and female peaks separately and find the highest peak in this "area" histogram and use that histogram to determine the y-coordinate value of the male or female peak.

[0195] 1. By using only events with event peaks close to the calculated "x" value, i.e. events with event peaks within 0.3% of 65536 (event peak range), the y coordinates of the male and female peak points where the events are most frequent can be found. The event peak range using the event area value to find the peak y coordinate is:

[0196] Peak value x-(65536*peak width×0.01) to

[0197] Peak value x + (65536*peak width×0.01).

[0198] 2. The event region values ​​for events close to the male peak and events close to the female peak are placed into separate arrays. Both arrays of "region" values ​​are sorted in ascending order.

[0199] 3. Create a histogram of the area values ​​for the male and female peaks respectively.

[0200] i. Use the following steps to calculate rather than predefine the bin widths for this histogram:

[0201] Calculate the number that is equal to 75% of the number of peak area values ​​(numValuesIn75pct);

[0202] Calculate the number that is equal to 99% of the number of peak area values ​​(largestRetainedIndex);

[0203] Ignore area values ​​whose values ​​are greater than the 99th percentile of all area values, as determined by the maximum retained index;

[0204] Determine the area value at the 99th percentile (maximum retained area);

[0205] Use numValuesIn75pct(lowest75pctMax) to find the area value of the 75th percentile of the marked area value;

[0206] Find the minimum area value in the sorted area value array from the events near the peak (lowest75pctMin);

[0207] Calculate the range of area values ​​in the lowest 75% of area values ​​(areaRange; = lowest75pctMax – minimum75pctMin);

[0208] · Calculate the bar width as areaRange / 512 (areaBinWidth);

[0209] Calculate the number of bins as (largestRetainedArea / areaBinWidth)+1(numberOfAreaBins).

[0210] If the number of area values ​​to be divided into bins is less than 6 times the number of bins just calculated, then numberOfAreaBins = (largestRetainedIndex / 6) + 1; if (numberOfAreaBins<4), numberOfAreaBins = 4; areaBinWidth = maximumRetainedArea / (numberOfAreaBins-1)

[0211] ii. After creating a histogram of area values ​​based on the calculated number of bins and bin widths, the histogram values ​​are smoothed using a 15-bin range centered around each smoothed bin value.

[0212] iii. Find the maximum value in the area histogram.

[0213] iv. Determine the area range and submit the least squares curve fitting algorithm.

[0214] Calculate the number of bars with an area that is 20% of the maximum area (20% of the "peak" area)

[0215] Moving leftward from the histogram's largest peak to the smaller region of bins, find the first bin with a count less than 20% of the largest bin count. This marks the leftmost bin submitted to the curve fitting algorithm.

[0216] Again, move rightward from the maximum value of the histogram to find the rightmost bar submitted to the curve fitting algorithm

[0217] v. Submit the area histogram to a least squares polynomial curve fitting algorithm and a Newton-Raphson algorithm to obtain more accurate area values ​​for the y-coordinate of each peak in the same units as the event area

[0218] If an error is generated while finding the y coordinate value, do not fall back to the highest point. The data is too sparse and even a smoothed histogram is too inaccurate to rely on the highest "area" bin as the peak. If a good "y" cannot be generated from the least squares polynomial curve fit, generate an error and abort the gate calculation

[0219] g. Calculate the female and male peak intensity coordinate values ​​and bin intensity values ​​adjusted for electronic device drift that has occurred since the event data set was processed from the most recent thirty event data sets. If fewer than 30 previous values ​​are available, then only the available values ​​are averaged.

[0220] 1. Create an event dataset counter variable and initialize it to zero.

[0221] 2. Each time a new set of event data is processed, the counter is incremented. The counter is incremented regardless of whether an error is flagged when processing the data and whether a character change in the event data is detected in step 8.

[0222] 3. When saving a value that will be averaged with similar values ​​in a subsequent data set for future use, save the value of the data set counter associated with the data set that produced the value.

[0223] 4. Calculate the average drift rate of variable X over the period of time spanned by processing multiple data sets by subtracting the earliest available saved value of X from the latest value of X, obtaining DELTA_X, and subtracting the value of the data set counter associated with the earliest value from the value of the counter associated with the most recently stored value, obtaining DELTA_T; then calculate DELTA_X / DELTA_T as the average drift.

[0224] 5. Calculate the total drift adjustment for a particular recorded instance of a variable by multiplying the average drift rate for the variable obtained from step d by the number of time periods that have passed since that value was recorded (i.e., the number of event data sets that have been presented). The product of the multiplications represents the drift adjustment. Add the drift adjustment to the previously recorded value to obtain the drift-corrected value.

[0225] 6. Calculate the peak-to-valley intensity distance by finding the difference between the bin intensity and the peak intensity of the peak where the gate was drawn.

[0226] 7. Generate an "average peak-to-valley intensity distance" by averaging the peak-to-valley intensity distance from the previous step and the peak-to-valley intensity distances from the most recent 30 event data sets. If fewer than 30 previous values ​​are available, only the available values ​​are averaged.

[0227] h. Calculate the slot x coordinate from the bars. Calculate the slot y coordinate so that the slot point is on a line drawn through the male and female peaks (x, y). Having determined the (x, y) coordinates of the male and female peaks, and determined the "x" coordinate of the slot bin between them, the "y" coordinate of the slot bin point is calculated to be on a line through the male and female peak points. This line is the axis of the gate.

[0228] i. Using the male peak, female peak coordinates and valley coordinates, calculate the gate vertex. The gate is oriented along the axis drawn between the male and female coordinates.

[0229] 1. Get the InvocementPCT value from the parameter set. Divide it by 100.0 and then subtract the result from 1.0 to get the encroachment multiplier used in the subsequent steps.

[0230] 2. Generate an encroachment adjustment by multiplying the average peak-to-valley intensity distance by the encroachment multiplier.

[0231] 3. Calculate the encroachment strength. The two vertices lie on a vertical line drawn through the encroachment strength.

[0232] i. If the gate is created relative to the male peak, the erosion intensity is the average male peak intensity plus the erosion adjustment. For the male peak, the erosion intensity is the intensity coordinate of the gate vertex or vertices that are the farthest to the right of the intensity axis in the (intensity, area) plane.

[0233] ii. If the gate is created relative to the female peak, the encroachment intensity is the average female peak intensity minus the encroachment adjustment. For the female peak, the encroachment intensity is the intensity coordinate of the leftmost gate vertex or vertices on the intensity axis in the (intensity, area) plane.

[0234] 4. Retrieve the peakFalloffPct value from the parameter set. This is a decimal value in the range 0.0 to 100.0. Currently 65.0 is used, but eventually some other value may be determined to be optimal. Divide peakFalloffPct by 100.0 and then multiply it by the height of the highest bar in the intensity histogram for the associated peak (step 10). The associated peak is either the male intensity peak or the female intensity peak, depending on which peak the gate is calculated relative to. In this process description, the product of the multiplication will be referred to as "falloffValue".

[0235] 5. Calculate the attenuation intensity: Starting from the highest intensity histogram bar in the correlation peak, move away from the bin in step 23 until the height of the currently examined bar is less than or equal to the attenuation value. Calculate the intensity value corresponding to the center of the bar as [bar index + 0.5] x [bar width].

[0236] 6. Calculate the average decay intensity by averaging the decay intensity of the previous step and the decay intensity of the previous thirty event data sets. If fewer than 30 previous values ​​are available, only average the available values.

[0237] 7. Take the average decay intensity as the initial value of the intensity coordinate in the plane of the two gate vertices (intensity, area), which are the rightmost in the gate drawn relative to the female peak, or the leftmost in the gate drawn relative to the male peak.

[0238] 8. Calculate the extension distance by retrieving the "octExtension" value from the parameter set and multiplying by the average peak-to-valley intensity distance. octExtension is a decimal value greater than 0.0. Use the default value of 0.6, but change it if desired.

[0239] i. If the gate is drawn relative to the convex peak, subtract the extension distance from the calculated initial average decay intensity value. Use the result as the final value of the intensity coordinate in the (intensity, area) plane for the two leftmost gate vertices in the gate drawn relative to the convex peak.

[0240] ii. If the gate is drawn relative to the parent peak, add the extension distance to the initial value. Use this result as the final value for the intensity coordinate in the (intensity, area) plane for the two rightmost gate vertices in the gate drawn relative to the parent peak.

[0241] 9. Calculate the mean invasion intensity and mean female peak intensity for gates drawn relative to the female peak, or the mean male peak intensity for gates drawn relative to the male peak.

[0242] 10. Retrieve the "horizadjmultipler1" value from the parameter set and multiply it by the average peak-to-valley intensity distance. Add the product to the average of the encroachment intensity values ​​calculated in the previous step.

[0243] i. For the gate relative to the female peak, the results are used as the intensity coordinates of the two vertices on the vertical line to the right of the encroachment intensity in the (intensity, area) plane.

[0244] ii. For gates relative to the male peak, the line is to the left of the eroded intensity. A value of 0.0 is currently used for HorizaDjMultipleR1, but this may change.

[0245] 11. Calculate the average of the mean decay intensity and the mean female peak intensity for gates drawn relative to the female peak, or the mean male peak intensity for gates drawn relative to the male peak.

[0246] 12. Retrieve the "horizadjmultipler2" value from the parameter set. Multiply it by the average peak-to-valley intensity distance. Add the product to the value calculated in the previous step.

[0247] i. For the gate relative to the female peak, the result is used as the intensity coordinate in the (intensity, area) plane of the two vertices on the vertical line to the right of the line calculated in step 10 and to the left of the line drawn through the decaying intensity.

[0248] ii. For the gating circuit relative to the male peak, the line is to the left of the line calculated in step 10 and to the right of the line plotted through the decaying intensities. A value of 0.0 is currently used for HorizaDjMultipler2 but this may change.

[0249] 13. Calculate the peak-to-valley area distance as the difference between the peak-to-valley "area" coordinates of the drawn gate and the peak-to-valley peak area coordinates.

[0250] 14. Calculate the average peak-to-valley area distance by averaging the peak-to-valley area distance from the previous step and the peak-to-valley area distances from the most recent 30 event data sets. If fewer than 30 previous values ​​are available, only the available values ​​are averaged.

[0251] 15. Calculate the intersection (intensity, area) where the axis from step 10 intersects the vertical line, and find the upper and lower vertices of the axis that fall on this line.

[0252] i. Calculate the "Area" coordinate for the vertex on the vertical line and above the axis by retrieving the value of "AreaAboverMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 14, then adding the result to the area coordinate of the intersection point.

[0253] ii. Calculate the "area" coordinates of the vertices below the axis in step 10 by retrieving the value of "areaBelowMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 14, then subtracting the result from the area coordinates of the intersection points.

[0254] 16. Position two vertices on a vertical line drawn through the encroachment point, one above and one below the axis line. Calculate the erosion intensity at the point where the vertical line intersects the axis line; in this process description, this point will be called the "encroachment point."

[0255] i. Calculate the “Area” coordinate above the vertex and axis line with the encroachment strength by retrieving the value of AboutHeightPctLeft from the parameter set, dividing it by 100.0, and then multiplying it by the length of the perpendicular line segment from the vertex to the axis line calculated in step 15; then add that value to the “Area” coordinate of the encroached point.

[0256] ii. Starting at the point where the vertical line at the erosion intensity intersects the axis line, calculate the vertex on the vertical line at the erosion intensity but below the axis line by retrieving the value of BelowEightPctLeft from the parameter set, dividing by 100.0, and multiplying by the length of the vertical line segment from the vertex to the axis line calculated in step 15; then subtract that value from the area value of the encroached point. This is the "area" coordinate of the point.

[0257] 17. Locate the two vertices that lie on the vertical line whose intensity was calculated in step 12. Calculate the intersection point (intensity, area) where the axis and the vertical line intersect in step 12.

[0258] i. Calculate the "Area" coordinate of the vertex on the axis line by retrieving the value of "AreaAboverMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 14, then adding the result to the area coordinate of the intersection point.

[0259] ii. Calculate the "area" coordinates of the vertices on the line below the axis from step 12 by retrieving the value of "areaBelowMultiplier" from the parameter set and multiplying it by the average peak-to-valley area distance from step 14, and subtracting the result from the area coordinates of the intersection points calculated for the vertices on the axis line.

[0260] 18. Position two vertices on a vertical line drawn through the falloff strength, one above the axis and one below the axis. Calculate the falloff point as the point where the vertical line at the falloff strength intersects the axis.

[0261] i. Calculate the Area coordinate of the vertex with decaying strength that is above the axis line by retrieving the value of AboutHeightPcTright from the parameter set, dividing it by 100.0, and multiplying it by the length of the perpendicular line segment from the vertex to the axis line calculated in step 17; then add that value to the Area coordinate of the decaying point.

[0262] ii. Starting from a falloff point with a falloff strength that is above the axis, calculate the area coordinate of the point whose vertex is on the vertical line at the falloff strength but below the axis by retrieving the value of BelowEightPCTright from the parameter set, dividing by 100.0, and multiplying by the length of the vertical line segment from the vertex to the axis calculated in step 17; then subtract that value from the area value of the falloff point. This is the "area" coordinate of the point below the axis.

[0263] Sex Identification of Bovine Sperm Cells Using a Program-Optimized Gating System

[0264] The laser ablation semen sexing cytometer was used to implement the optimized sorting strategy described in the above example. Semen was collected from a 2-year-old Holstein bull and prepared according to standard procedures. Briefly, the collected semen was evaluated for motility and the presence of abnormalities to find that it met the processing criteria. Sperm cell concentration was then determined and the sperm was stained with Hoechst 33342 for 45 minutes at 34°C before sexing. The sexing instrument includes a flow controller to deliver the stained sperm cells and sheath fluid to the microfluidic chip. Sperm cell detection uses TM Mode-locked 355 laser The fluorescence was detected by an avalanche photodiode. The laser was Q-switched using the ExplorerXP laser. Ablation of unwanted sperm cells carrying the Y chromosome. The dyed perm cells were treated for approximately 9 hours.

[0265] The discrimination of sperm cells carrying X- and Y-chromosomes was performed using an optimized gating strategy according to an embodiment of the present invention, as detailed in the above example. After classification, the collected sample produced approximately 51.46 million sperm cells. During operation, the values ​​of each gate parameter were monitored and the operation of the gate system was observed. During a single sample run, extracted portions of the intensity value histogram generated from events at two time points ten minutes apart were compared. The curves were fitted to peaks and valleys; mathematical analysis of the fitted curves yielded more accurate peak and valley values ​​displayed above each figure. The exact peak and valley values ​​were converted into actual intensity values ​​by multiplying the bar size (128). Thus, for example, the male peak intensity in the histogram on the left is 223.06×128.0=28551.68.

[0266] Comparison of the peak intensity values ​​at the two different time points shows that the female peak intensity changes significantly from 29625 to 29386 (239 intensity units, about 0.81%) in 10 minutes. The bin also moves from 29101 to 28832 (269 intensity units). To maintain the correct relationship to the female peak position, the gate must move by a similar amount.

[0267] Since the gate is built around the female peak, similar histograms are programmatically constructed for the "Area" values ​​of all events with intensities close to the intensity of the female peak. These histograms show only one peak. The analyzed histograms of the "Area" values ​​show that the female peak area has a peak value of 4972333 at the first time point and a peak value of 4952019 at the second time point.

[0268] All intensity coordinates of the gate shift leftward (lower) between two time points within ten minutes. The intensity coordinate of the gate apex at the first time point averaged 29735.28. The corresponding average for the second time point (10 minutes later) was 29465.78. This corresponds to an average downward shift of 269.5 intensity units. This compares favorably to the female peak intensity shift of 239 intensity units and the trough intensity shift of 269 intensity units. The gate shifted leftward to maintain its relationship with the female peak.

[0269] The area coordinates of the female peak also shifted downwards over the 10 min period. The area coordinates of all gates at the apex shifted downwards in a similar manner.

[0270] The discrepancy between the exact amount of displacement of the instantaneous peak coordinates and the gate apex coordinates is mainly due to the use of the trough-female peak distance in the calculation of the apex. Since the trough and female peak are independent statistical entities, the apex coordinates may shift slightly differently than the female peak between two time points, despite the fact that the average shifts of the apex and peak coordinates over time are equal.

[0271] Those skilled in the art will appreciate that various changes and substitutions of equivalent technologies may be made without departing from the scope of the new technologies disclosed in this application. In addition, many modifications may be made to adapt specific situations or materials to the teachings of the new technologies without departing from their scope. Therefore, the new technologies are not limited to the specific technologies disclosed, but include all technologies that fall within the scope of the appended claims.

Claims

1. A method for identifying a subpopulation in a particle mixture, the method being used to identify a subpopulation in a particle mixture, comprising: measuring a characteristic of a first group of particles in a mixture of particles having two or more subpopulations, each subpopulation being distinguishable by the characteristic; Measured characteristic values ​​that distinguish each subpopulation in a mixture of particles; Based on the distinction, a cutoff range is assigned to each subpopulation; measuring a characteristic of a second group of particles in the particle mixture after assigning the critical value range; After measuring the characteristic of the second group of particles in the particle mixture, distinguishing the measured characteristic value of each subpopulation in the second group of particles; updating the critical range of values ​​based on distinguishing the measured characteristic values ​​of each subpopulation in the second group of particles; and assigning subpopulation classifications to particles in the particle mixture according to the updated critical value range; wherein the particles are cells; The characteristics include size, shape, permeability, integrity, content, viability, mortality, metabolic state, surface markers, DNA content, number of events, presence of an X chromosome or a Y chromosome, or a combination thereof.

2. The method according to claim 1, further comprising: Based on the updated threshold range, the particle is classified according to the classification assigned to the particle.

3. A method for identifying multiple groups of particles in a mixture of particles, each particle subgroup having characteristics different from other particle subgroups, comprising: providing a mixture of particles; detecting a characteristic of a first group of particles in the particle mixture to generate a signal indicative of the characteristic of the first group of particles; generating event data for a signal indicative of a characteristic of the first set of particles; Compile event data for the first set of particles; Generate a distribution of event data for the first set of particles; calculating at least one discriminant value defining a subpopulation of particles based on a distribution of event data of the first group of particles; Classifying particles based on discriminant values ​​that define particle subpopulations; wherein the particles are cells; The characteristics include size, shape, permeability, integrity, content, viability, mortality, metabolic state, surface markers, DNA content, number of events, presence of an X chromosome or a Y chromosome, or a combination thereof; The event data includes fluorescence intensity, total fluorescence, or a combination thereof.

4. The method according to claim 3, characterized in that: Also included is classifying the device to include or exclude events based on the classification indication.

5. The method according to claim 3, characterized in that: Calculating at least one discriminant value includes defining a critical value associated with at least the detected feature.

6. The method according to claim 5, characterized in that The calculating of at least one discriminant value further comprises: identifying at least one critical value of a first feature from a distribution of event data for a first set of particles; obtaining a value of a second feature associated with at least one critical value of the first feature; Producing a distribution of second eigenvalues; identifying at least one critical value of the second feature from the distribution of values ​​of the second feature; co-locating the critical values ​​of the first and second features in the feature space; and At least one discriminant value is determined from critical values ​​that are co-located in the feature space.

7. The method according to claim 6, characterized in that The calculating of the at least one discriminant value further includes adjusting the thresholds of the first feature and the second feature to include event data from a second set of particles.

8. The method according to claim 6, characterized in that The critical values ​​for the first feature and the second feature include one or more peak feature values ​​for the first particle subpopulation and one or more peak feature values ​​for the second particle subpopulation.

9. The method according to claim 8, characterized in that The peak characteristic values ​​include peak fluorescence intensity and peak total fluorescence.

10. The method according to claim 8, characterized in that The critical values ​​of the first characteristic and the second characteristic also include valley values ​​between the peak characteristic values ​​of the first and second particle subpopulations.

11. The method according to claim 3, characterized in that The at least one discriminant value comprises a gate.

12. The method of claim 1 or 3, further comprising destroying or disabling a subpopulation of particles assigned the classification.

13. The method of claim 1 or 3, wherein the particle is a sperm cell, wherein the sperm cell is characterized by the presence of an X chromosome or a Y chromosome.

14. A sample of sexed sperm cells produced by the method of any one of claims 1-13.

15. A system for producing a sexed sperm sample, comprising: a sample delivery device for delivering the sperm cell mixture, the sample delivery device delivering a sample of the sperm cell mixture to a flow cytometer system; a flow cytometer for measuring characteristics of a first group of sperm cells in the sperm cell mixture, wherein the characteristics can distinguish different subpopulations; The control device is: The measured characteristic values ​​that distinguish each subpopulation in a mixture of sperm cells, Based on the distinction, a critical value range is assigned to each subgroup. updating the critical range of values ​​based on distinguishing the measured characteristic values ​​of subpopulations in the second group of sperm cells, and assigning subpopulation classifications to sperm cells in the mixture based on the updated cutoff range; and A sorting device, wherein the sorting device performs a sorting operation on sperm cells in a sperm cell mixture according to instructions received from the control device; the system is used to implement the method described in any one of claims 1-2 and 12-14.

16. The system according to claim 15, characterized in that The sorting device includes a laser killing device.

17. The system according to claim 15 or 16, wherein the control device can be further configured to: identifying at least one critical value for the first feature from the distribution; obtaining a value of a second feature associated with at least one critical value of the first feature; Producing a distribution of second eigenvalues; identifying at least one critical value of the second feature from the distribution of the second feature values; jointly locating a critical value of the first feature and the second feature in the feature space; as well as At least one discriminant value is determined from critical values ​​of a first feature and a second feature co-located in the feature space.