A method for quantitatively evaluating the spreading of a rhizomatous plant

Abstract

The present application relates to the technical field of computer-aided modeling and intelligent agricultural phenotype analysis, and particularly relates to a rhizomatous plant expansion quantification evaluation method, comprising: obtaining rhizome imaging data and environmental variable data of the same rhizomatous plant to be evaluated at at least two sampling time points; determining a set of expansion front nodes from a rhizome skeleton graph and constructing a rhizome digital twin model; establishing a set of front actions oriented to the set of expansion front nodes; constructing a consistency objective function; determining a confidence index based on the quality of observation data and the stability of registration, and adaptively determining a weight coefficient according to the confidence index; performing iterative calibration on a parameter vector based on the consistency objective function, and performing constraint projection on the updated parameter vector in each iteration to meet a constraint set; calculating an expansion index based on a predicted state trajectory and outputting an evaluation result. The present application realizes stable and quantitative evaluation of rhizome expansion.

Description

Technical Field

[0001] This invention relates to the field of computer-aided modeling and intelligent agricultural phenotypic analysis technology, and in particular to a method for quantitative evaluation of the extensibility of rhizomatous plants. Background Technology

[0002] The underground organs of rhizomatous plants (such as plants with distinct rhizomes or stolons) typically exhibit structural evolutionary behaviors during growth, including extension along the leading edge, local branching, turning in different directions, and termination under specific conditions. This structural evolution is reflected not only in geometrical changes (e.g., changes in length, radius, and spatial distribution) but also in topological changes (e.g., the addition of branch nodes, branch edges, and changes in connection relationships). Therefore, in agricultural breeding, cultivation management, and phenotypic research, there is an urgent need to quantitatively assess the "expansion" of rhizomatous plants to enable comparative analysis of different samples under different environmental conditions.

[0003] In existing technologies, the digital structural characterization of root systems is typically established by imaging underground organs (e.g., X-ray CT) and reconstructing their three-dimensional structures (e.g., reference 1: four-dimensional root structure modeled from a digital twin using X-ray computed tomography). Based on this, structural traits such as total length, number of branches, diameter, and volume are extracted for analysis. For example, published literature has proposed a four-dimensional (time-varying three-dimensional) digital twin modeling method for root systems based on X-ray CT. This method achieves temporal analysis and quantitative research of root structures by performing skeletal characterization, registration, and trait extraction on the structures at different times.

[0004] However, based on the practical needs of quantitative assessment of the spreadability of rhizomatous plants, the aforementioned existing technologies may still face the following objective problems in engineering applications:

[0005] First, the imaging quality, degree of occlusion or missing parts, density of reconstructed points, and skeleton connectivity vary at different sampling times, which affects the stability and comparability of the evaluation results when cross-time comparisons are made based on a single geometric metric or single trait statistics.

[0006] Secondly, the rootstock expansion process includes discrete growth behaviors such as "extension, branching, turning, and termination". The corresponding structural evolution involves both geometric morphological changes and topological structural changes. Although existing technologies can obtain skeletal characterization and extract traits, they lack a unified quantification and calibration mechanism for "the coupling consistency between topological changes and geometric changes during structural evolution". This makes it difficult to simultaneously take into account geometric matching, topological matching, and phenotypic matching in prediction or comparative analysis.

[0007] Third, in the presence of noise, fractures and missing parts, the evaluation process often needs to meet structural rationality constraints such as bifurcation angle, minimum spacing of skeleton nodes, and diameter conservation. Existing technologies lack a mechanism to use the above constraints as feasible domain conditions for consistency calibration during parameter updates or model calibration, which affects the repeatability of scalability evaluation.

[0008] Therefore, the main technical problems that need to be solved by existing technologies can be summarized as follows: how to establish a calibrable digital twin evaluation process based on the root and stem skeleton structure under the condition that the quality of imaging data fluctuates at multiple time points and the evolution of root and stem structure includes both geometric and topological changes. By uniformly constraining geometric consistency, topological consistency and phenotypic consistency and completing calibration within the feasible domain that meets the structural rationality constraints, a stable quantitative evaluation and cross-time comparable analysis of the extensibility of root and stem plants can be achieved. Summary of the Invention

[0009] To overcome the above-mentioned technical defects, the present invention aims to provide a quantitative assessment method for the extensibility of rhizomatous plants. The present invention uses digital twin prediction of rhizomatous skeleton diagram driven by frontier action sequence, and performs closed-loop calibration under confidence adaptive weighting and constrained projection calibration through a unified objective function composed of geometric consistency, topological editing action consistency and phenotypic consistency, thereby achieving a stable quantitative assessment of rhizomatous extensibility.

[0010] This invention discloses a method for quantitatively evaluating the extensibility of rhizomatous plants, characterized by being executed by electronic devices and comprising the following steps:

[0011] S1, acquire root and stem imaging data and environmental variable data of the same rhizomatous plant to be evaluated at at least two sampling times; perform three-dimensional reconstruction on the root and stem imaging data at each sampling time to obtain the observation point set. and Rhizome skeleton maps were generated based on the observation point set. and ,in, The sampling time number, A set of skeleton nodes. Create a skeleton edge set; configure node attributes for the skeleton nodes of the rhizome skeleton graph and edge attributes for the skeleton edges. Node attributes include at least node coordinates and node radius attributes, and edge attributes include at least edge length attributes.

[0012] S2, from the rhizome skeleton diagram Determine the set of extended frontier nodes And construct a digital twin model of the rhizome, defining the model state as Define the environment input as Define the model parameters as a parameter vector. ;

[0013] S3, Establishing a set of nodes for expanding the frontier The set of leading-edge actions, which includes at least: extension actions, bifurcation actions, turning actions, and termination actions; based on environmental input. With parameter vector To expand the set of frontier nodes Generate a leading action sequence and apply the leading action sequence to the rhizome skeleton diagram. Obtain the predicted skeleton map ;

[0014] S4, Construct the consistency objective function Consistency objective function Includes: a set of observation points Compared with the predicted skeleton diagram Geometric consistency term obtained from the generated prediction point set 1. Rhizome skeleton diagram With predicted skeleton map Topology consistency terms obtained from differences in topology editing actions And the phenotypic consistency term derived from the difference between the observed phenotypic vector and the predicted phenotypic vector. And satisfy:

[0015]

[0016] in, , These are the weighting coefficients;

[0017] S5. The confidence index is determined based on the quality of the observation data and the stability of registration, and the weighting coefficients are adjusted according to the confidence index. , Perform adaptive determination;

[0018] S6, based on the consensus objective function For parameter vectors Perform iterative calibration, and in each iteration, update the parameter vector. Perform constraint projection to satisfy the constraint set Constraint Set It includes at least: minimum spacing constraints for skeleton nodes, bifurcation angle constraints, and diameter conservation constraints; and is based on the calibrated parameter vector. Generate predicted state trajectories;

[0019] S7, Calculating Scalability Indicators Based on Predicted State Trajectory It also outputs the evaluation results, including the scalability index. At least with the extended frontier node set The displacement increment, coverage area increment, and volume increment are related.

[0020] Preferably, in step S1, the rhizome skeleton map is obtained by performing rhizome region segmentation and skeleton extraction on the observation point set.

[0021] Preferably, the generation of the rhizome skeleton map further includes performing skeleton pruning on the skeleton extraction results. Skeleton pruning includes at least: deleting isolated short branches and performing connection completion on broken skeleton segments to form skeleton connected components.

[0022] Preferably, in step S2, the set of leading edge nodes is expanded. By screening root and stem skeleton diagrams The terminal node with a degree of 1 is obtained, and the root node corresponding to the starting point of the rootstock is removed from the set of nodes at the expanding front edge. Excluded from the list.

[0023] Preferably, the end node that intersects with the preset imaging space boundary is marked as the boundary front node, and boundary constraints are applied to the boundary front node to limit its action type or action parameters when generating the front action sequence in step S3.

[0024] Preferably, in step S3, the set of leading edge actions includes: an extension action that updates the position of the leading edge node along the extension direction, a forking action that generates at least one sub-branch of the leading edge node, a turning action that updates the extension direction of the leading edge node, and a termination action that stops the growth of the leading edge node.

[0025] Preferably, the leading-edge action sequence is represented in the form of a triple of "leading-edge node - action type - action parameter", wherein: the action parameters of the extension action include at least the extension step size and extension direction, the action parameters of the forking action include at least the number of sub-branches and the set of sub-branch directions, the action parameters of the turning action include at least the direction rotation angle, and the action parameters of the terminating action include at least the termination marker and the termination time index.

[0026] Preferably, the forking action is defined by a forking rule base, which includes: the conditions for allowing forking at the leading edge node, the range of possible forking child nodes, and the range of possible initial directions for the child branches; and at least a portion of the rule parameters in the forking rule base are defined by a parameter vector. Provided.

[0027] Preferably, in step S4, the topology consistency item By drawing the rhizome skeleton diagram With predicted skeleton map These are converted into topology editing action sequences and aligned to obtain the results.

[0028] Preferably, the topology editing actions include at least: node insertion, node deletion, edge insertion, edge deletion, and node attribute replacement; node attribute replacement includes at least replacing the node radius attribute and the edge length attribute; and the cost of the topology editing action is determined based on the action type cost and the attribute difference cost, and the topology editing action sequence corresponding to the minimum total cost is used to determine the topology consistency item. .

[0029] Preferably, in step S5, the confidence index includes at least: geometric confidence, topological confidence, and registration confidence; wherein, the geometric confidence is determined by the observation point set. The point density and missing rate were determined, and the topological confidence was obtained from the rhizome skeleton map. The number of broken segments and the number of connected components are determined, and the registration confidence is determined by the registration residuals of adjacent sampling times; the registration residuals are the observation point sets of adjacent sampling times. and The residuals obtained from the registration were performed; and the weighting coefficients were adjusted according to the confidence index. , Adaptive determination is performed.

[0030] Preferably, the weighting coefficients , The adaptive determination adopts discrete rules: when the number of broken segments exceeds the preset threshold, the lower limit of the weight coefficient corresponding to the topological consistency term is set to be no less than the preset lower limit value; when the point missing rate exceeds the preset threshold, the upper limit of the weight coefficient corresponding to the geometric consistency term is set to be no greater than the preset upper limit value; and the weight coefficient is normalized under the condition of satisfying the upper and lower limit constraints.

[0031] Preferably, in step S6, the diameter conservation constraint limits that the sum of the squares of the radii of the multiple sub-branches obtained by bifurcation at the same parent node does not exceed the square of the parent node's radius; the bifurcation angle constraint limits that the angle between the newly added skeleton edge and its parent skeleton edge falls within a preset angle range; and the minimum spacing constraint limits that the distance between any two different skeleton nodes is not less than a preset minimum spacing.

[0032] Preferably, the parameter vector Including continuous and discrete regular parameters, the constraint projection includes at least: performing boundary trimming on continuous parameters, performing nearest-neighbor replacement within the feasible value set on discrete regular parameters, and correcting bifurcation action parameters that violate diameter conservation constraints; and regenerating the frontier action sequence after completing the constraint projection to update the predicted skeleton graph. .

[0033] Preferably, the predicted state trajectory is generated using a rolling time-domain method. When root and stem imaging data at newly added sampling times are obtained, the parameter vector is adjusted. Perform backtracking calibration on the values ​​at least two most recent sampling times, and then convert the backtracked parameter vector... Used for generating leading-edge action sequences in the subsequent time domain; and the rhizome skeleton diagram and calibrated parameter vectors are then used. Consistency objective function Iteration records, scalability metrics Write to the storage unit to generate a structured evaluation report.

[0034] Compared with existing technologies, the above technical solution has the following advantages:

[0035] 1. Achieve quantifiable and traceable assessment of the "expansion" of rhizomatous plants: By acquiring rhizome imaging data and environmental variable data at at least two sampling times, a consistent set of observation points and rhizome skeleton map across time periods are formed, and the expansion index, its components, and traceability information are output, so that the expansion assessment is transformed from a qualitative description into a calculable and verifiable quantitative result.

[0036] 2. Maintaining evaluation stability under imaging defects, noise, and structural fracture conditions: By using a multi-threshold judgment and minimum cost matching mechanism for skeleton trimming, fracture determination, and fracture completion, the impact of skeleton fractures and misconnections caused by imaging holes and local defects on structural characterization is reduced, thereby enhancing cross-time structural consistency and evaluation stability.

[0037] 3. Achieve reliable calibration of the twin model through the three consistency objectives of "geometric-topological-phenotypic": Construct a consistency objective function that includes geometric consistency terms, topological consistency terms and phenotypic consistency terms, and integrate constraint space fit, structural connectivity and phenotypic statistical consistency during the iterative calibration process to form an optimizable closed-loop consistency relationship between the established rhizome digital twin model and the observed state.

[0038] 4. Improving robustness and applicability through adaptive weighting based on confidence index: Applying upper and lower bound constraints and normalizing the weights of consistency items based on quality markers such as missing rate, number of broken segments, number of connected components, and registration residuals, so that the objective function can automatically adjust its focus under different data quality and structural complexity conditions, reducing the bias caused by a single consistency item.

[0039] 5. Ensuring the rationality of biological morphology and suppressing unreal growth through constrained projection closed loop: After parameter update, constraints such as minimum spacing of skeleton nodes, bifurcation angle and diameter conservation are projected and corrected. Based on the projected parameters, the leading-edge action sequence and predicted skeleton diagram are regenerated, reducing structural extrapolation that does not conform to morphological rules and improving the reliability of predicted trajectory and expansion index.

[0040] 6. Clear phenotypic calculation criteria facilitate cross-sample and cross-condition comparisons: Phenotypic quantities such as total length, number of branches, number of extension front nodes, coverage area, spatial extension radius and volume are calculated from the skeleton diagram using a unified criterion, and the coverage area and volume estimation methods are given with standard mathematical formulas, improving the consistency of comparisons between different batches, different imaging scales and different culture conditions.

[0041] 7. The results are highly interpretable and easy to verify and apply in engineering: By recording traceability information such as sampling time, environmental window statistics, registration matrix and registration residual, break completion records, weight coefficient values, objective function iteration sequence and constraint projection trigger records, the scalability evaluation results can be located to the specific data source and processing link, which is convenient for verification, adjustment and systematic deployment.

[0042] 8. Facilitates integration with existing processes and supports different imaging and acquisition conditions: The solution supports multiple imaging sources and environmental acquisition configurations at the data acquisition level, and provides examples of implementable parameters and experimental configuration methods in tabular form, enabling it to be migrated and implemented under different equipment conditions and application scenarios, reducing integration and deployment costs. Attached Figure Description

[0043] Figure 1 This is a schematic diagram illustrating the steps of a quantitative assessment method for the extensibility of rhizomatous plants according to the present invention;

[0044] Figure 2 This is a diagram of the rhizome skeleton, an expanded set of frontier nodes, and frontier actions.

[0045] Figure 3 This is a schematic diagram of the closed loop of weight adaptation and constraint projection;

[0046] Figure 4 This is a schematic diagram of the convergence curve of the consistency objective function;

[0047] Figure 5 This is a schematic diagram illustrating the relationship between the missing rate and the stability of the expansion index.

[0048] Figure 6 This is a schematic diagram illustrating the consistency of the scalability index scatter plot. Detailed Implementation

[0049] The advantages of the present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments.

[0050] Exemplary embodiments will now be described in detail, examples of which are illustrated in the accompanying drawings. When the following description relates to the drawings, unless otherwise indicated, the same numerals in different drawings denote the same or similar elements. The embodiments described in the following exemplary embodiments do not represent all embodiments consistent with this disclosure. Rather, they are merely examples of apparatuses and methods consistent with some aspects of this disclosure as detailed in the appended claims.

[0051] The terminology used in this disclosure is for the purpose of describing particular embodiments only and is not intended to be limiting of the disclosure. The singular forms “a,” “the,” and “the” as used in this disclosure and the appended claims are also intended to include the plural forms unless the context clearly indicates otherwise. It should also be understood that the term “and / or” as used herein refers to and includes any and all possible combinations of one or more of the associated listed items.

[0052] It should be understood that although the terms first, second, third, etc., may be used in this disclosure to describe various information, such information should not be limited to these terms. These terms are used only to distinguish information of the same type from one another. For example, without departing from the scope of this disclosure, first information may also be referred to as second information, and similarly, second information may also be referred to as first information. Depending on the context, the word "if" as used herein may be interpreted as "when," "when," or "in response to determination."

[0053] In the description of this invention, it should be understood that the terms "longitudinal", "lateral", "up", "down", "front", "rear", "left", "right", "vertical", "horizontal", "top", "bottom", "inner", "outer", etc., indicate the orientation or positional relationship based on the orientation or positional relationship shown in the accompanying drawings. They are only for the convenience of describing this invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, or be constructed and operated in a specific orientation. Therefore, they should not be construed as limitations on this invention.

[0054] In the description of this invention, unless otherwise specified and limited, it should be noted that the terms "installation", "connection" and "linking" should be interpreted broadly. For example, they can refer to mechanical or electrical connections, or internal connections between two components. They can be direct connections or indirect connections through an intermediate medium. Those skilled in the art can understand the specific meaning of the above terms according to the specific circumstances.

[0055] In the following description, suffixes such as "module," "part," or "unit" used to denote elements are used only for the convenience of the description of the invention and have no specific meaning in themselves. Therefore, "module" and "part" can be used interchangeably.

[0056] This embodiment provides a method for quantitatively assessing the extensibility of rhizomatous plants, executed by an electronic device, including the following steps: S1, acquiring rhizome imaging data and environmental variable data of the same rhizomatous plant to be assessed at at least two sampling times; performing three-dimensional reconstruction on the rhizome imaging data at each sampling time to obtain an observation point set. and Rhizome skeleton maps were generated based on the observation point set. and ,in, The sampling time number, A set of skeleton nodes. S1: Create a skeleton edge set; configure node attributes for the skeleton nodes and edge attributes for the skeleton edges of the rhizome skeleton graph. Node attributes must include at least node coordinates and node radius attributes, and edge attributes must include at least edge length attributes; S2: From the rhizome skeleton graph... Determine the set of extended frontier nodes And construct a digital twin model of the rhizome, defining the model state as Define the environment input as Define the model parameters as a parameter vector. S3, Establish a set of nodes for extended frontiers. The set of leading-edge actions, which includes at least: extension actions, bifurcation actions, turning actions, and termination actions; based on environmental input. With parameter vector To expand the set of frontier nodes Generate a leading action sequence and apply the leading action sequence to the rhizome skeleton diagram. Obtain the predicted skeleton map S4, Construct the consensus objective function Consistency objective function Includes: a set of observation points Compared with the predicted skeleton diagram Geometric consistency term obtained from the generated prediction point set 1. Rhizome skeleton diagram With predicted skeleton map Topology consistency terms obtained from differences in topology editing actions And the phenotypic consistency term derived from the difference between the observed phenotypic vector and the predicted phenotypic vector. And satisfy: ,in, , S5 represents the weighting coefficients; based on the quality of the observed data and the stability of the registration, the confidence index is determined, and the weighting coefficients are adjusted according to the confidence index. , Perform adaptive determination; S6, based on the consistency objective function For parameter vectors Perform iterative calibration, and in each iteration, update the parameter vector. Perform constraint projection to satisfy the constraint set Constraint Set It includes at least: minimum spacing constraints for skeleton nodes, bifurcation angle constraints, and diameter conservation constraints; and is based on the calibrated parameter vector. S7: Generate predicted state trajectory; S8: Calculate scalability index based on predicted state trajectory. It also outputs the evaluation results, including the scalability index. At least with the extended frontier node set The displacement increment, coverage area increment, and volume increment are related.

[0057] This specific embodiment provides a method for quantitatively assessing the extensibility of rhizomatous plants, executed by an electronic device. The electronic device includes a processor, a memory, an imaging data interface, an environmental data interface, and a result output interface. The processor calls program instructions from the memory to process rhizome imaging data and environmental variable data, outputting extensibility assessment results, extensibility indicators, and their traceability information. For ease of explanation, this specific embodiment uses X-ray CT imaging and multi-sensor environmental acquisition under root box cultivation conditions as an example implementation environment, but this does not limit the imaging method or environmental variable acquisition method; in other embodiments, rhizome imaging data can also be obtained by industrial CT, structured light point cloud, multi-view imaging reconstruction, or ultrasound data, and environmental variable data can also be obtained by other measuring devices.

[0058] To avoid confusion in technical terms and for ease of understanding, this embodiment provides symbolic identifiers for some quantization objects: observation point set. With observation point set Indicates: Rhizome skeleton diagram. With rhizome skeleton diagram Representation; Extended front node set is used to represent the extended front node set. Indicates; predicting skeleton map using predicting skeleton map Indicates; the consistency objective function is used. Indicates; geometric consistency terms are used. Indicates; topology consistency terms are used. Indicates; Phenotypic consistency terms are used. Indicate; weighting coefficients are used. Weighting coefficients Weighting coefficients Representation; constraint sets are used to represent constraint sets. Indicates; Expansion indicators are used. express.

[0059] like Figure 1 As shown, the method in this embodiment includes the following steps within one evaluation cycle: data acquisition and time window alignment; obtaining an observation point set through 3D reconstruction; calibration marker identification and cross-time registration; generating a root and stem skeleton map from the observation point set and performing pruning and break completion; determining the set of extended frontier nodes and constructing a root and stem digital twin model; generating a frontier action sequence and obtaining a predicted skeleton map; constructing a consistency objective function containing geometric consistency, topological consistency, and phenotypic consistency terms, adaptively determining weight coefficients based on confidence indices, and performing iterative calibration and constrained projection closure; calculating phenotypic vectors and extensibility indices and outputting the results; and performing comparative experiments and statistical tests to demonstrate the effectiveness and stability when verification is required.

[0060] In this embodiment, a root box is used as a growth carrier for underground rhizomes. The interior of the root box is filled with a matrix and maintained at a controllable moisture content. At least three spatial calibration markers are fixedly installed on the outside of the root box. The spatial calibration markers appear as high-density small spherical structures or structures with significant contrast in the imaging data, and are used for coarse registration and attitude consistency verification across sampling times. The root box is fixed in an imaging fixture, which is equipped with an attitude positioning structure to ensure that the root box is installed in a consistent attitude at different sampling times, thereby reducing rotational and translational deviations in cross-time registration.

[0061] The environmental variable acquisition system includes at least one soil moisture sensor, at least one soil temperature sensor, and at least one soil conductivity sensor. The sensors sample at fixed intervals via a data logger and upload the data with timestamps to an electronic device. The electronic device stores the sampling record of each variable as a data row of "timestamp-variable name-variable value-quality marker" and uses the environmental variable window statistics corresponding to each imaging time as environmental input.

[0062] The electronic device establishes a unique identifier for each sample and creates a sampling record for each sampling. The sampling record includes at least: sampling time, imaging data path, environmental window statistics, calibration marker identification results, registration matrix, observation point set and rootstock skeleton map pointers, and intermediate result pointers required for subsequent consistency and scalability index calculations. This organizational method ensures that scalability indexes are traceable to specific sampling times, specific imaging data, and specific environmental windows. See Tables 1 and 2 for details.

[0063] Table 1: Examples of Imaging Equipment and Scanning Parameters

[0064]

[0065] Table 1 illustrates how key scanning parameters of the imaging equipment are set when acquiring root and stem imaging data, ensuring that the volumetric data output has spatial resolution and signal-to-noise ratio suitable for subsequent 3D reconstruction, segmentation, and point set generation. In Table 1, "tube voltage and tube current" define X-ray energy and irradiation intensity, which together affect the grayscale contrast and noise level of the volumetric data; "voxel size" provides the spatial sampling interval for the 3D volumetric data, directly determining the dimensional basis for subsequent point cloud sampling step size, minimum identifiable structural scale of the skeleton, and threshold range; "rotation step size and number of projection frames" characterize the angular coverage and sampling density of projection acquisition, affecting the reconstruction artifact level; "reconstruction algorithm and output format" describe the volumetric data generation method and file readability, ensuring that the electronic device can generate 3D volumetric data based on metadata stacking or projection reconstruction and enter the subsequent process. The example values ​​given in Table 1 visually demonstrate that this embodiment adopts an implementable and reproducible scanning configuration in the "data acquisition stage," providing a scale basis for subsequent calculations such as "point density, missing rate, and threshold range."

[0066] Example Table 2 of Environment Variable Collection Configuration

[0067]

[0068] Table 2 illustrates the source, sampling rhythm, and window statistics configuration of environmental variable data collected in conjunction with root and stem imaging data. In Table 2, "Variable Name" lists the types of environmental variables involved in this embodiment; "Sensor Location / Number" indicates that multiple sampling points can be used for the same variable to reduce sporadic noise and enhance representativeness; "Sampling Period" determines the temporal resolution of the environmental variable data, thereby determining the theoretical number of sampling points within a certain window length; "Window Length" specifies the time span for extracting the environmental window near the sampling time, facilitating the calculation of statistics such as mean, standard deviation, and missing rate; "Dimensions" defines the physical units of the variable to ensure interpretability of subsequent thresholds and comparisons; and "Example of Quality Marker" indicates the reliability evaluation field of the variable window, such as missing rate reflecting data gaps, and Std reflecting short-term fluctuations or noise. Through the configuration in Table 2, it is clear that this embodiment enters environmental variables into the electronic device in the form of a "time-stamped continuous sampling sequence," and provides a basis for subsequent action parameter adjustment, weight adaptation, and traceability recording in the form of "window statistics + quality markers."

[0069] In this embodiment, to ensure comparability across time periods, records from two sampling times of the same rhizomatous plant to be evaluated constitute one evaluation cycle. The two sampling times can specifically be sampling at the same time on the same day for two or more days at intervals. In this embodiment, 09:00 on the 7th day of cultivation is used as the sampling time. The sampling time was 09:00 on the 10th day of cultivation. The electronic device sends a scanning task to the imaging device at the sampling time. The scanning task includes at least: sample identifier, scanning parameter template, output path, and scanning sequence number. After the root box is installed in the imaging fixture, the imaging device performs a flat-field acquisition for subsequent correction, followed by sample scanning to acquire projection data or slice data. After scanning, the imaging file is output. The electronic device reads the output directory through the imaging data interface and performs verification on the continuity and integrity of the file sequence. For example, it verifies whether the slice sequence number is consecutive, whether the number of slices meets expectations, and whether the voxel size and scanning parameters in the imaging metadata are consistent with the task template. After the verification is successful, the imaging file pointer is written to the sampling record.

[0070] Environmental variable data contain transient noise and occasional missing data around the sampling time. Therefore, this embodiment extracts a fixed-length environmental window for statistical analysis at each sampling time. Specifically, this can be implemented as follows: for each sampling time... Extract window data from 30 minutes before the sampling time to the sampling time; for the sampling time Similarly, data from 30 minutes prior to the sampling time to the sampling time are extracted. The electronic device calculates the mean, standard deviation, and missing rate for each environmental variable within the window, and uses the missing rate and standard deviation as quality markers. Taking soil moisture content as an example, when the sampling period is 5 minutes, the 30-minute window theoretically contains 6 points. If one point is missing, the missing rate is 1 / 6. When the missing rate exceeds a preset threshold (e.g., 20%), the electronic device marks the variable window as low confidence and reduces the sensitivity of the action parameters driven by this variable or increases the lower limit of the topology consistency weight in subsequent weight adaptation.

[0071] To ensure consistent attitude across sampling times, this embodiment employs both mechanical positioning and calibration mark verification. Mechanical positioning is achieved using an imaging fixture, ensuring a unique attitude after root box insertion. Calibration mark verification involves identifying the center coordinates of calibration marks in the volumetric data. If the relative position of the identified calibration mark differs from the reference attitude by more than a threshold, the electronic device marks this sampling as an attitude anomaly and prompts for reinstallation or the addition of stability in subsequent registration. See Table 3 for details.

[0072] Example Table 3 of Environmental Window Statistics Fields at Two Sampling Times

[0073]

[0074] Table 3 illustrates how environmental window statistics are presented under the same window rules for two sampling times, and how these statistics are written into the sampling record and used for subsequent processing in this embodiment. In Table 3, "Sampling Time, Window Start and End" defines the time boundaries for window extraction; "Mean" represents the representative level of environmental variables within the window; "Std" represents the intensity of short-term fluctuations within the window; and "Missing Rate" represents the proportion of missing sampling points within the window. The example in Table 3 demonstrates that even with consistent sampling times (e.g., both at 09:00), window statistics from different dates may still differ. Electronic devices can use this to generate environmental input and synchronously record quality markers. When the missing rate exceeds a threshold or Std increases abnormally, it can trigger upper and lower limit constraints on subsequent weighting coefficients or reduce sensitivity to environmental drivers, thereby avoiding amplified effects of poor environmental data quality on the modeling results.

[0075] At the sampling time respectively With sampling time After obtaining the root and stem imaging data, the electronic device performs three-dimensional reconstruction at each sampling time to generate a set of observation points. This process can be implemented as a continuous link of the following sub-steps.

[0076] First, when the input is projection data, the electronic device performs flat-field correction and dark-field correction to eliminate system bias, performs ring artifact removal to suppress ring artifacts, and performs filtering to reduce high-frequency noise. Then, filtered backprojection is used to reconstruct the 3D volume data. When the input is a slice sequence, the electronic device stacks the slices according to metadata to form 3D volume data and normalizes and denoises the voxel intensity, making the intensity distribution at different sampling times comparable. Denoising can be achieved using methods such as 3D median filtering or nonlocal mean filtering to preserve the root and stem boundaries.

[0077] Secondly, the electronic device performs root and stem region segmentation on the 3D volumetric data. Root and stem region segmentation can be specifically implemented as follows: a candidate voxel set is obtained using threshold segmentation, where the threshold can be determined by an automatic thresholding method or given by calibration samples; subsequently, morphological closing operations are performed on the candidate voxel set to fill small holes, and connected component screening is performed to retain larger connected components while eliminating small-volume noisy connected components; when there are high-brightness spots in the matrix particles, the spatial range of the candidate region can be further limited using spatial calibration markers and root box boundaries, thereby reducing the probability of missegmenting non-root and stem structures as root and stem structures. The electronic device records the segmented voxel volume, the number of connected components, and the proportion of hole area as segmentation quality indicators.

[0078] Next, the electronic device performs surface extraction on the segmented root and stem voxel regions. Surface extraction can be specifically implemented by generating a triangular mesh using an isosurface extraction algorithm, and then simplifying the triangular mesh to reduce redundant faces. Simplification can be achieved by edge folding while maintaining the overall shape error below a threshold. Subsequently, the electronic device performs uniform sampling on the surface of the triangular mesh to generate an observation point set. Uniform sampling can be specifically implemented as Poisson disk sampling, ensuring that the minimum spacing between sampling points is not less than a preset value, thereby guaranteeing uniform point distribution; alternatively, it can be implemented as weighted random sampling based on face area and downsampling of locally overly dense areas. Each point in the observation point set contains at least three-dimensional coordinates, and if necessary, a normal vector and intensity value. The electronic device calculates the number of points, average neighborhood distance, and point density for the observation point set, and estimates the missing rate, which can be estimated through the proportion of mesh hole area or through point cloud coverage. The aforementioned point density and missing rate are used for subsequent geometric confidence calculation and weight adaptation.

[0079] In a feasible recording example, the set of observation points at sampling time k The number of points is approximately 185,000, and the sampling time is... observation point set The number of points is approximately 172,000, and the sampling time is... The missing rate is estimated to be 20%. The electronic device writes this missing rate into the sampling record and uses it as the trigger condition for subsequent weight constraints.

[0080] In this embodiment, the spatial calibration marker can be a high-density spherical marker or a high-contrast structural marker. The electronic device's identification of the calibration marker in 3D volume data can be specifically implemented as follows: First, a threshold screening is performed on high-intensity voxels in the volume data to obtain a candidate voxel set. Then, connected component extraction is performed on the candidate voxel set to obtain several candidate connected components. For spherical markers, the electronic device calculates the number of voxels, the ratio of principal inertia along the three axes, and the sphericity index for each candidate connected component, and selects connected components that meet the criteria of "the number of voxels falling within the marker volume range and the sphericity being higher than the threshold" as marker candidates. For each marker candidate, the electronic device calculates its center coordinates. The center coordinates can be obtained by calculating the centroid of the connected component voxel coordinates, or by fitting a sphere within the candidate region and taking the center coordinates of the sphere, thereby improving the accuracy of the center coordinates. The electronic device writes the center coordinates of each marker as the calibration marker identification result into the sampling record.

[0081] To avoid misidentification, the electronic device can further check whether the relative distance between the calibration marks is consistent with the preset mark layout. When the difference between the center distance of any two marks and the reference distance exceeds the threshold, the electronic device marks the identification as unreliable and can use stable estimation or prompt re-imaging in the next registration step.

[0082] At the time of sampling With sampling time After determining the center coordinates of the calibration markers, the electronic equipment solves for the coarse registration rigid body transformation across time intervals. Specifically, this can be implemented by: [The text abruptly ends here, likely due to an incomplete sentence or missing information.] The set of three marker center coordinates is denoted as the reference marker set, and the sampling time is... The set of coordinates of the centers of the three markers is denoted as the target marker set. Rigid body registration in the least squares sense is used to solve for the rotation matrix and translation vector, ensuring that the target marker set is as close as possible to the reference marker set after the rigid body transformation. The solution can be obtained using a closed-form method based on singular value decomposition, forming a coarse registration matrix. The electronic device applies this coarse registration matrix to the sampling time. The set of observation points or volume data, making it consistent with the sampling time. Enter the same coordinate system.

[0083] Coarse registration primarily eliminates large-scale attitude differences, but millimeter-level errors may still exist. To improve the reliability of the consistency term calculation, this embodiment further performs fine registration. Fine registration can be specifically implemented as follows: for sampling times... Observation point set Sampling time after coarse registration Observation point set The iterative nearest neighbor algorithm is executed. To improve stability, the electronic device can first perform voxel downsampling on the point cloud to reduce the number of points and calculate the normal vector of the point cloud to support the point-to-surface error measurement. In each iteration, the electronic device finds the nearest neighbor of the reference point cloud for each point in the target point cloud, calculates the point-to-surface error and solves for a small rigid body increment, and updates the registration matrix. Iteration stops when the error decrease is less than a threshold or the number of iterations reaches the upper limit. To suppress the influence of outliers, the electronic device can remove matching point pairs with a nearest neighbor distance greater than a threshold, or use a truncation error function to reduce the weight of outliers.

[0084] After fine registration, the electronic device calculates the registration residual. Specifically, the registration residual can be calculated by taking the root mean square of the Euclidean distance or point-to-surface distance for all retained matching point pairs and using this root mean square as the registration residual. The registration residual is used to form the registration confidence score; the smaller the registration residual, the higher the registration confidence score. The electronic device writes the registration matrix and registration residual into the sampling record and uses the registration residual as one of the trigger conditions during subsequent weight adaptation. For example, when the registration residual is greater than a threshold, the lower bound of the topology consistency term weight is increased to reduce the sensitivity of the geometric consistency term to errors.

[0085] The generation of a root and stem skeleton map from an electronic device using an observation point set can be specifically implemented as follows: First, the observation point set is pixelated to generate an occupying voxel set, the voxel size of which can be the same as or slightly larger than the imaging voxel; then, centerline extraction is performed on the occupying voxel set, which can be obtained by three-dimensional refinement to obtain centerline voxels; then, the centerline voxels are converted into a graph structure, the adjacency relationship between the centerline voxels is converted into skeleton edges, and the bifurcation points or degree changes are defined as skeleton nodes, thereby forming a skeleton node set and a skeleton edge set.

[0086] After skeleton extraction, short branches may exist due to noise. The electronic device can specifically remove short branches by calculating the path length from each terminal node to the nearest branch node, marking terminal branches with path lengths less than a threshold as short branches and deleting them. The threshold can be determined based on the voxel size and the minimum identifiable root segment length of the target; for example, the threshold can be in the range of 5 mm to 15 mm. To avoid deleting genuine fine root branches, the electronic device can also combine node radius attributes for judgment. If the average radius attribute of a short branch is greater than the threshold and its directional continuity is good, the branch is retained; otherwise, it is deleted.

[0087] Skeleton fractures typically manifest as the presence of multiple connected components in the skeleton graph, or the existence of suspected fracture endpoints within the same connected component. The electronic device can determine fracture endpoints by: selecting end nodes with a degree of 1 as candidate endpoints in the skeleton graph; calculating the local direction vector for each candidate endpoint, which can be obtained by taking several nodes inward along the skeleton edge from that endpoint and fitting a straight line; statistically analyzing the voxel density or point cloud density of the connected component containing the candidate endpoint; if the point cloud density near the candidate endpoint significantly decreases and is adjacent to matrix pores, then the candidate endpoint is marked as a fracture endpoint. The electronic device counts the number of fracture endpoint pairs to obtain the number of fracture segments, and also counts the number of connected components as input for topological confidence.

[0088] The key to fracture completion lies in "when to connect, which to connect, how to connect, and what shape to connect." This embodiment uses a multi-threshold judgment and minimum cost matching method to reduce false connections. Specifically, it can be implemented as follows: for all fracture endpoints, form candidate pairs, and calculate the distance between endpoints, directional consistency, radius continuity, and boundary crossing risk.

[0089] Distance threshold determination: When the distance between endpoints is greater than the maximum completion distance threshold, no completion is performed. The maximum completion distance threshold can be set according to the voxel size and imaging defect scale, for example, within the range of 3 mm to 20 mm, and can be appropriately increased as the defect rate increases.

[0090] Directional consistency determination: Calculate the angle between the local direction vectors of the two endpoints. If the angle is greater than the maximum angle threshold, do not complete the calculation. The maximum angle threshold can be taken in the range of 20 degrees to 60 degrees.

[0091] Radius continuity determination: Calculate the difference in node radius attributes between the two endpoints. When the radius difference exceeds the threshold, reduce the completion priority or disable completion. The radius difference threshold can be taken in the range of 0.2 mm to 0.6 mm.

[0092] Boundary crossing judgment: If the completed line crosses the preset imaging space boundary or crosses a high-density non-root area, the completion is prohibited.

[0093] To take into account the above factors, this embodiment defines a completion cost for each candidate pairing. The completion cost is a weighted sum of the distance, direction, and radius terms. The electronic device selects the non-conflicting pairing with the minimum total cost from all candidate pairs that meet the hard threshold and performs completion. Here, "non-conflicting" means that the same broken endpoint cannot be paired with multiple endpoints at the same time.

[0094] After determining the endpoint pairings, the electronic device can generate completed skeleton edges in one or a combination of the following ways: First, use linear interpolation to insert several intermediate nodes between the two endpoints and connect them to form completed edges; second, use spline fitting, using the endpoint positions and local direction vectors as boundary conditions to fit a smooth curve and insert intermediate nodes along the curve; third, perform a shortest path search on the endpoints in the voxel occupancy set to obtain pathways, and map the pathways to skeleton edges. To avoid unreasonable sharp turns caused by completed edges, the electronic device can perform smoothing processing on local nodes after completion, such as performing moving average or curvature constraint smoothing on node coordinates, and recalculate the edge length attribute.

[0095] After completion, the electronic device updates the number of connected components in the skeleton graph and records the difference in the number of broken segments before and after completion as a completion effect indicator. Simultaneously, the completion trigger record is written into the traceability information. See Table 4 for details.

[0096] Example Table 4: Key Thresholds and Recommendation Range

[0097]

[0098] Table 4 illustrates the semantics, dimensions, and recommended ranges of the thresholds used in this embodiment when dealing with key steps such as "pruning, break completion, registration, and outlier suppression." The "short branch deletion threshold, maximum completion distance threshold, maximum included angle threshold, and radius difference threshold" in Table 4 are used to constrain the triggering conditions and connectivity determination of skeleton pruning and break completion, thereby reducing noisy branches and erroneous completions. The "registration residual threshold" is used to evaluate the reliability of cross-time registration and provide triggering conditions for subsequent weight adaptation. The "nearest neighbor truncation threshold" is used to suppress the influence of outliers on distance statistics during geometric consistency term calculation. Table 4 emphasizes the "feasible threshold caliber," whose recommended range matches the voxel scale, point cloud density, and missing rate level, enabling those skilled in the art to select reasonable thresholds and reproduce the skeleton completion and consistency calculation process at a given imaging resolution.

[0099] The node radius attribute represents the equivalent radius of the root section at that node. Specifically, the electronic device calculates the node radius attribute as follows: Centered on the coordinates of the skeleton node, a neighborhood point cloud or voxel set is extracted in a plane perpendicular to the local direction vector. The set of distances from the neighborhood points to the node center is recorded as the candidate radius set. The median or truncated mean is used as the equivalent radius to reduce the influence of outliers. The neighborhood truncation radius can be in the range of 1 mm to 5 mm, and the truncation thickness can be in the range of one to multiple voxel thicknesses. The electronic device stores the obtained equivalent radius as the node radius attribute.

[0100] The edge length attribute represents the length of a skeleton edge. Specifically, the electronic device calculates the edge length attribute by summing the node sequence contained in a skeleton edge according to the coordinate distances of adjacent nodes, and then using this path length as the edge length attribute. For edges generated by break completion, the electronic device also sums the interpolated node sequence to obtain the edge length attribute, ensuring consistent length calculation.

[0101] Electronic devices at sampling time The set of nodes for the expansion front is determined on the rhizome skeleton diagram. Specifically, this can be implemented as follows: terminal nodes with a degree of 1 are selected as candidate terminal nodes; root nodes are identified and excluded from the candidate terminal nodes. Root node identification can be implemented as follows: determining the upper edge direction in the root box coordinate system, selecting terminal nodes near the upper edge that belong to the path with the largest diameter as root nodes, or manually inputting the root node positions during the initial modeling of the sample and mapping the root node positions along the registration matrix at subsequent times. The electronic device uses the remaining terminal nodes as the set of nodes for the expansion front.

[0102] For boundary constraints, after determining the set of extended leading edge nodes, the electronic device calculates the minimum distance from each extended leading edge node to the preset imaging space boundary. When the minimum distance is less than the safe distance threshold, the extended leading edge node is marked as a boundary leading edge node, and the minimum distance is recorded. When generating subsequent action sequences, the electronic device performs directional projection on the extension direction of the boundary leading edge nodes or performs clipping on the extension step size, so that the newly added nodes are still within the boundary and maintain a safe distance, thereby reducing evaluation distortion caused by exceeding the boundary.

[0103] like Figure 2 As shown, the electronic device establishes a set of front-end actions, including at least extension actions, bifurcation actions, turning actions, and termination actions. The electronic device generates an action sequence for each extension front node, which can be specifically implemented as follows: first, based on environmental input and model parameter vectors, determine the action type probability or action triggering condition; then, based on the action type, determine the action parameters. Environmental inputs can be used to adjust the extension step size and bifurcation probability; for example, reducing the extension step size or increasing the termination probability when water content decreases. The model parameter vectors provide adjustable coefficients and thresholds, such as the extension step size coefficient, turning strength coefficient, and bifurcation trigger threshold.

[0104] To ensure that forking behavior is interpretable and constrainable, this embodiment establishes a forking rule base. The forking rule base includes at least: a threshold for the node radius attribute that allows for forking, a threshold for the minimum trunk length that allows for forking, a set of possible sub-branch numbers, a range of initial direction angles for sub-branchs, and sub-branch radius allocation rules. When generating a forking action, the electronic device first checks whether the expansion front node meets the allowed forking conditions. If it does, it takes a value from the set of possible sub-branch numbers and generates a set of sub-branch directions within the included angle range. The sub-branch radius allocation rules can be consistent with the diameter conservation constraint, ensuring that the sum of the squares of the sub-branch radii does not exceed the square of the parent radius during initial generation. If subsequent iterations cause a violation of the constraint, it is corrected by constraint projection.

[0105] The electronic device writes the action sequence into the tracing information in the format of "Extended Front Node Identifier - Action Type - Action Parameter - Time Index," and applies the actions to the rootstock skeleton graph to obtain the predicted skeleton graph. When applying an extension action, new nodes and edges are added; when applying a forking action, multiple child nodes and edges are added, and the extended front node set is updated; when applying a turning action, the orientation state of the extended front node is updated, affecting subsequent extension directions; when applying a termination action, the extended front node is removed from the extended front node set or set to an inactive state. For boundary front nodes, the electronic device performs boundary constraint correction before applying action parameters, and the correction record is also written into the tracing information.

[0106] The key to topology consistency lies in "topology editing action alignment". This embodiment treats the stem-and-skeleton graph as an approximate tree structure with the root node as the root, since the stem-and-skeleton graph is in most cases a tree-like or approximate tree structure. The electronic device can convert the skeleton graph into a root tree by performing a breadth-first traversal starting from the root node, determining the parent node of each node and forming parent-child relationships, thus obtaining the root tree structure; when local loops caused by noise exist, edges with larger lengths or radii can be retained first during traversal, and loop edges can be deleted, thus converting the structure into a tree.

[0107] To ensure stable tree editing alignment, this embodiment performs a consistent sorting of the child node set for each node. Specifically, the sorting can be implemented as follows: using the parent edge direction as a reference, calculate the azimuth angle of each child edge direction in the reference plane and sort by azimuth angle. If the azimuth angles are similar, then sort by child edge length, thus obtaining an ordered list of child nodes. This sorting ensures that the same type of branching structure has a more consistent child node order at different sampling times, facilitating dynamic planning alignment.

[0108] After sorting the root tree and child nodes, the electronic device performs tree edit distance calculation on the two trees. The tree edit distance calculation can be specifically implemented as a minimum-cost alignment based on dynamic programming: allowed edit actions include node insertion, node deletion, edge insertion, edge deletion, and node attribute replacement. Node insertion and deletion correspond to structural additions and deletions, edge insertion and deletion correspond to connection relationship additions and deletions, and node attribute replacement corresponds to adjusting the difference in node radius or edge length attributes. The electronic device defines an action type cost for each type of action and an attribute difference cost for attribute replacement. The attribute difference cost can be obtained by multiplying the absolute difference of the attributes by a scaling factor. The difference in node radius attribute is converted to a dimensionless cost using "millimeters multiplied by a factor," and the difference in edge length attribute is also converted to a dimensionless cost using "millimeters multiplied by a factor." The electronic device calculates the minimum cost for each pair of node subtrees in dynamic programming and accumulates it from bottom to top to obtain the minimum total cost of the entire tree. This minimum total cost is the topology consistency term. The electronic device can also output the corresponding editing action sequence as traceability information, making the source of topological differences explainable. For example, it can output "Insert a node and an edge at a certain fork, and replace the node radius attribute and the edge length attribute".

[0109] In a sample implementation, the cost of node insertion can be 10, the cost of edge insertion can be 8, and the cost of node and edge deletion can be the same as or slightly lower than the insertion cost. The cost of node attribute replacement consists of the base cost of 2 plus the attribute difference cost, where the difference in node radius attribute is multiplied by 20, and the difference in edge length attribute is multiplied by 5. The above cost values ​​are used to demonstrate an implementable cost system, and can be adjusted according to the data scale in actual engineering.

[0110] Phenotypic vectors are calculated item by item: the phenotypic consistency term requires both the observed phenotypic vector and the predicted phenotypic vector. The calculation of phenotypic vectors from the rhizome skeleton diagram in electronic devices can be specifically implemented by including at least one or more of the following phenotypic quantities, and the corresponding calculation methods are given.

[0111] Total length calculation: The total length is the sum of the length attributes of all skeleton edges in the root-stem skeleton graph. The electronic device traverses the skeleton edge set and sums the edge length attributes of each edge to obtain the total length. If the skeleton is a tree structure, the total length can also be regarded as the sum of the lengths of all parent-child connection paths.

[0112] Branch Count Calculation: The branch count can be specifically implemented as the number of forked nodes or the number of forked edges. A forked node can be defined as a node with a degree of not less than 3. The electronic device counts the number of nodes that meet the condition as the branch count; or the sum of the number of sub-branches of each forked node is counted as the branch edge count. To maintain consistency, this embodiment uses the same metric for comparison.

[0113] Calculation of the number of extended front nodes: The number of extended front nodes is the number of elements in the set of extended front nodes. The electronic device directly counts these elements after determining the set of extended front nodes.

[0114] Coverage area calculation: Coverage area measures the spatial extent of the rootstock's spread on a projection plane. The electronic device can calculate the coverage area as follows: First, determine the projection plane, which can be the root box side plate plane, the root box bottom plate plane, or the principal plane obtained from principal component analysis of the skeleton; then, project the skeleton node coordinates or the coordinates of the observation point set onto this projection plane to obtain a two-dimensional projection point set; finally, establish a mesh on the two-dimensional plane, with the mesh side length being the mesh step size. The grid step size can be taken from 0.5 mm to 5 mm and adjusted according to the voxel size. For each grid cell, it is determined whether it is covered by a projection point. This determination can be implemented as follows: if there is at least one projection point within the grid cell, it is recorded as an occupied cell; or if the number of projection points within the grid cell exceeds a threshold, it is recorded as an occupied cell. Finally, the covered area is obtained by multiplying the number of occupied cells by the area of ​​the grid cell. Its standard expression can be written as:

[0115]

[0116] The number of occupied units is the number of occupied units. The grid step size is the grid step size. To reduce the underestimation of coverage area caused by holes, electronic devices can perform a two-dimensional closing operation on the occupied grid to fill small holes, and only fill holes with an area smaller than a threshold to avoid filling real blank areas.

[0117] Calculation of spatial expansion radius:

[0118] The spatial expansion radius measures the maximum outward expansion of a root-and-stem structure from the root node. Specifically, calculating the spatial expansion radius in electronic devices involves calculating the Euclidean distance from all skeleton nodes to the root node or the distance along the skeleton path, and taking the maximum value as the spatial expansion radius. Using the distance along the skeleton path makes it more sensitive to structural bending; using the Euclidean distance provides a more intuitive view of the overall spatial extent.

[0119] Volume calculation:

[0120] The volume is calculated from the skeleton edges and node radii. Electronic devices can treat each skeleton edge as a segment with a variable cross-section, approximating the volume as the sum of the volumes of segmented frustums or cylinders. To avoid using special symbols, this embodiment uses the constant pi. In electronic devices, the radius attributes of the nodes at both ends of a skeleton edge are respectively denoted as the endpoint radius attributes. With endpoint radius attribute The side length attribute is denoted as the side length attribute. Then the volume of the frustum corresponding to that side can be approximated as:

[0121]

[0122] The electronic device iterates through all skeleton edges and accumulates the sum of each edge. The total volume is obtained. If the skeleton sides are short or the radius does not vary much, the average radius can also be used. Approximate volume of a cylinder:

[0123]

[0124] Electronic devices maintain a consistent volume calculation standard within the same system, and the standard used is written into the traceability information.

[0125] Growth rate calculation:

[0126] The growth rate can be obtained by dividing the difference between the total length or volume at two sampling times by the time interval. For example, the growth rate of the total length is... ,in This represents the time difference between two sampling points. Electronic devices can include growth rate as one of the phenotypic quantities in the phenotypic consistency term calculation. As shown in Table 5:

[0127] Table 5 Summary of Phenotypic Vector Items and Calculation Methods

[0128]

[0129] Table 5 serves as an example, illustrating which phenotypic quantities are included in the "phenotypic vector" in this embodiment, and the calculation object, calculation method, and unit for each phenotypic quantity. The "total length, number of branches, and number of expansion front nodes" in Table 5 are primarily obtained from the node and edge attributes of the rhizome skeleton graph; "coverage area" emphasizes the calculation method for the area occupied after projection and meshing; "spatial expansion radius" reflects the maximum expansion scale outward from the root node; "total volume" is used to combine the node radius attribute with the edge length attribute to obtain a volume estimate; and "growth rate" is used to obtain the dynamic change quantity through time-separation. Table 5 clearly and verifiably fixes the mapping relationship "from skeleton to phenotypic," avoiding inconsistencies caused by different algorithms or units for the same phenotypic quantity in different sections; it also facilitates the selection of dimensions and the calculation of normalized differences in the phenotypic consistency item.

[0130] The electronic device calculates the geometric consistency term by converting the predicted skeleton graph into a predicted point set. Specifically, it samples along each skeleton edge at a fixed step size, ranging from 0.5 mm to 3 mm and matching the voxel size. Simultaneously, it records the corresponding local radius attribute for each sampled point (obtained through linear interpolation of the radii at both ends of the edge), forming a predicted point set with radius attributes. The electronic device then calculates the average nearest neighbor distance from the observed point set to the predicted point set, and vice versa. Nearest neighbor queries can be implemented using a spatial index structure, such as constructing a kd-tree to accelerate the query. To reduce the impact of outliers, the electronic device can remove or truncate point pairs with nearest neighbor distances exceeding a threshold, and write the removal ratio into a quality flag. The geometric consistency term is the sum or weighted sum of the aforementioned bidirectional average distances; this embodiment uses the sum of bidirectional average distances to reduce the bias caused by unidirectional missing values.

[0131] The calculation of phenotypic consistency terms by electronic devices can be implemented as follows: First, select the dimensions of the phenotypic vector, such as three or more of the following: total length, total volume, number of branches, coverage area, and number of nodes in the extended front. For each dimension, determine a normalization benchmark, which can be the value of that dimension observed or a preset scale. Calculate the absolute difference for each dimension and divide it by the normalization benchmark to obtain the normalized difference. Sum the normalized differences according to preset weighting coefficients to obtain the phenotypic consistency term. If a phenotypic quantity is unusable due to its absence (e.g., the projection quality of the coverage area is too low), the electronic device can remove that phenotypic quantity from the vector and re-normalize the remaining weights, or set the weight of that phenotypic quantity to 0 and record the reason for its absence to maintain the computation's feasibility.

[0132] The electronic device can calculate confidence indices as follows: Geometric confidence is determined by the density and missing rate of the observation point set; high geometric confidence occurs when the density is high and the missing rate is low. Topological confidence is determined by the number of connected components and broken segments in the skeleton; high topological confidence occurs when the number of connected components and broken segments is low. Registration confidence is determined by the registration residual; high registration confidence occurs when the registration residual is less than a threshold. The electronic device maps these confidence indices to weighted constraint rules, which can use discrete thresholds for ease of implementation and review. For example, when the missing rate is greater than 10%, the upper limit of the geometric consistency term weight coefficient is limited to no more than 0.35; when the number of broken segments is greater than 3 or the number of connected components is greater than 2, the lower limit of the topological consistency term weight coefficient is limited to no less than 0.40; when the registration residual is greater than 1.5 mm, the upper limit of the geometric consistency term weight coefficient is further reduced and the lower limit of the topological consistency term weight coefficient is increased. After satisfying the upper and lower limit constraints, the electronic device selects a set of feasible weights and performs normalization processing so that the sum of the weight coefficients is 1, thereby maintaining the scale stability of the consistency objective function. As shown in Table 6:

[0133] Example Table 6 of Weight Adaptive Rules

[0134]

[0135] Table 6 illustrates how electronic devices apply upper and lower bound constraints to weight coefficients based on confidence-related conditions when observation data is missing, the skeleton is broken, or the registration residual is large. This ensures that the contribution ratio of each consistency term in the consistency objective function matches the data quality. The "missing data rate condition" in Table 6 limits the upper limit of the geometric consistency term weight coefficient to avoid the geometric distance statistics being misled by missing structures when observations are highly missing. The "number of broken segments / number of connected components condition" increases the lower limit of the topological consistency term weight coefficient to strengthen the constraints on structural connectivity and edit alignment. The "registration residual condition" further suppresses geometric terms and strengthens topological terms when cross-time alignment is unstable, thereby reducing the amplification effect of registration errors on the objective function. Table 6 demonstrates how the triggering logic of "weight adaptation" can be written as an implementable rule set, making the weight determination process interpretable and reproducible, and forming a consistent link with iterative calibration and constraint projection closed loop.

[0136] The consistency objective function is defined as:

[0137]

[0138] The iterative calibration of the model parameter vector by the electronic device can be implemented as follows: An feasible optimization strategy is selected, such as coordinate descent, random perturbation search, or gradient approximation update based on finite difference. Taking finite difference approximation as an example, the electronic device applies a small perturbation to each continuous parameter in the model parameter vector, calculates the difference in the consensus objective function before and after the perturbation, estimates the improvement direction in that parameter direction, and updates the parameter according to the step size. Discrete regular parameters can be updated by enumerating in the candidate set and selecting the value that minimizes the consensus objective function. After each round of updates, the electronic device records the value of the consensus objective function and terminates the process using the maximum number of iterations or the magnitude of the decrease in the objective function, for example, stopping when the difference in the consensus objective function between two adjacent rounds is less than a threshold or the number of iterations reaches the upper limit. Figure 4 The curve representation of this type of iterative process is given, showing that the consistency objective function decreases and stabilizes with iteration.

[0139] Constraint projection is used to ensure that the updated parameters and the structure generated from the parameters satisfy biological morphological rationality. The constraint set includes at least the minimum spacing constraint between skeleton nodes, the bifurcation angle constraint, and the diameter conservation constraint. The minimum spacing constraint between skeleton nodes can be achieved by checking whether the distance between any two adjacent new nodes is less than a threshold; if it is less than the threshold, the nodes are merged or their positions are reverted. The bifurcation angle constraint can be achieved by calculating the angle between the newly added edge and the parent edge and trimming it to within the threshold range. The diameter conservation constraint can be expressed by the following formula:

[0140]

[0141] The parent node radius attribute is the parent node radius attribute. The child node radius attribute is the child node radius attribute. The number of sub-branches is the number of sub-branches. When the diameter conservation constraint is violated, the electronic device can use scaling to correct the child node radius attribute, making the sum of squares of the child node radii equal to or less than the sum of squares of the parent node radius. The scaling factor can be obtained by multiplying the child node radii by the scaling factor and ensuring the sum of squares satisfies the constraint. After constraint projection, the electronic device must regenerate the leading-edge action sequence and update the predicted skeleton diagram based on the projected model parameter vector to ensure consistency between parameters, actions, and structure; this closed-loop process is as follows: Figure 3 As shown, the reason for each projection trigger and the correction amount are written into the traceability information.

[0142] After completing calibration and obtaining the predicted state trajectory, the electronic equipment calculates the expansion index. The expansion index can be specifically implemented as a comprehensive index related to displacement increment, coverage area increment, and volume increment. Displacement increment can be defined as the average or total displacement of the set of nodes at the expansion front before and after the predicted advance. The electronic equipment calculates the displacement for each node at the expansion front and averages or sums the displacement increment. Coverage area increment can be defined as the difference in coverage area between two time points, calculated using the aforementioned projected grid aperture. Volume increment can be defined as the difference in total volume between two time points, calculated using the aforementioned truncated cone accumulation plus aperture. The electronic equipment can obtain the expansion index using a linear combination method.

[0143]

[0144] Where the displacement increment is the displacement increment The increase in coverage area is the increase in area. The volume increment is the volume increment. The coefficients a, b, and c are non-negative and can be given by a preset parameter or by the model parameter vector, remaining non-negative in the constraint projection. When the electronic device outputs the scalability index, it also outputs traceability information. This traceability information includes at least: sampling time, registration residual, missing rate, number of broken segments, weight coefficient values, consistency objective function iteration sequence, constraint projection trigger record, and components of the scalability index. , , And its calculation caliber identifier.

[0145] To verify the stability and consistency of this embodiment under conditions of missing data and complex structures, a comparative experiment is further provided. The comparative experiment can be specifically implemented as follows: Rhizomatous plant samples of the same variety and under the same culture conditions are selected as experimental samples, with a sample size of 30 or more. Sample inclusion criteria may include: no obvious disease before sampling time, no significant displacement of the root box, complete imaging data file and continuous slice sequence, and identifiable calibration markers that meet distance consistency checks. Sample exclusion criteria may include: severe motion artifacts in the imaging leading to segmentation failure, calibration marker recognition failure or abnormal relative distance between markers, and registration residuals exceeding a preset threshold and failing to converge through fine registration. After including the samples, the electronic device can randomly number the samples to avoid human selection bias and perform the same acquisition process on each sample at both sampling times.

[0146] The comparison method can be a method that does not include the key closed-loop mechanism of this embodiment, such as a modeling process based solely on geometric registration and trait difference, without employing topological editing alignment, weight adaptation, and constrained projection regeneration closed loops. The electronic device calculates the geometric matching error, topological difference cost, scalability index stability, and scalability index consistency for both the method of this embodiment and the comparison method. The scalability index stability can be specifically implemented by repeatedly performing reconstruction and skeleton extraction three times on the same sample. The repetition method can be changing the segmentation threshold perturbation, changing the surface sampling random seed, or changing the skeleton refinement initialization, thereby obtaining three scalability index results, and calculating the coefficient of variation as the stability index. The scalability index consistency can be specifically implemented by manually labeling the expansion area increment of the sample as a reference quantity, calculating the correlation coefficient between the scalability index and the reference quantity, and expressing it as a scatter plot, as shown in the figure. Figure 6 As shown.

[0147] To verify robustness under missing rate conditions, the comparative experiment can be specifically implemented as a missing rate simulation: a certain proportion of points are randomly deleted from the observation point set to form a missing observation point set. The missing rate can be 0%, 10%, 20%, or 30%. For each missing rate condition, the stability of the extensibility index is repeatedly calculated, and the relationship curve between missing rate and stability is plotted. The curve is presented as follows: Figure 5 As shown in Table 7, statistical tests can be implemented as follows: a paired test is used for the stability index obtained from the same sample under both methods; a paired t-test is used when the index distribution is approximately normal, and a Wilcoxon signed-rank test is used when the distribution does not meet the normality requirement; the significance level can be set to 0.05, and the p-value is output to indicate the statistical significance of the difference. The electronic device writes the sample list of the comparative experiment, the reasons for inclusion and exclusion, the missing rate simulation parameters, the replication and reconstruction strategy, and the statistical test results into the experimental record for verification.

[0148] Example Table 7 of Comparative Experiment Design and Statistical Test Configuration

[0149]

[0150] Table 7 serves as an example, illustrating how experimental samples are selected, repeated reconstructions are set up, missing rate simulations are implemented, and statistical tests are conducted when verifying the stability and consistency of the method in this embodiment under different missing rate conditions. The "Sample Size, Inclusion Criteria, and Exclusion Criteria" in Table 7 clarify the source of experimental data and quality control principles, ensuring that the comparison results are not dominated by obvious artifacts or unusable samples. "Number of Repeated Reconstructions" defines the number of repetitions and variable factors for stability assessment, thereby calculating the coefficient of variation. "Missing Rate Simulation" involves randomly deleting observation points to construct missing conditions and generate curves. "Stability Indicators and Consistency Indicators" specify the statistical objects. "Statistical Tests" describes the significance assessment methods and significance levels for differences between comparative methods. Through Table 7, this embodiment provides a complete chain of data selection, experimental control, and statistical testing in an engineering-executable manner, making the method's effectiveness verification verifiable and reliable.

[0151] For ease of understanding, the following example is provided in this embodiment:

[0152] In the example record, the sampling time Sampling time: 09:00 on day 7 of cultivation. Sampling time: Day 10, 09:00 The observation point set has a missing rate of 20%, four broken segments, and a registration residual of 1.3 mm. The electronic equipment sets the upper limit of the geometric consistency term weight coefficient to 0.35 and the lower limit of the topological consistency term weight coefficient to 0.40 according to the rules. Feasible unnormalized weight combinations of 0.35, 0.45, and 0.30 are selected and then normalized to obtain the weight coefficients. Weighting coefficients Weighting coefficients .

[0153] The electronic device discretely samples the predicted skeleton map to obtain a predicted point set and calculates the average bidirectional nearest neighbor distance. For example, the nearest distances from six observation points to the predicted point set are 0.8, 1.0, 0.9, 1.2, 0.7, and 1.1 mm, yielding an average of 0.95 mm; the nearest distances from predicted points to the observation point set are 0.6, 0.8, 0.7, 0.9, 0.6, and 0.8 mm, yielding an average of 0.73 mm. (Geometric consistency term...) The value is 1.68. The electronic device performs tree editing alignment on the observed and predicted skeletons. An example optimal editing action sequence includes one node insertion, one edge insertion, and two attribute replacements. The topology consistency term is obtained according to the action cost system. For calculating phenotypic vector differences in electronic devices, an example is to take the total length, total volume, and number of branches, and then apply normalized weights to obtain the phenotypic consistency term. The electronic device calculates the consistency objective function based on this and uses iterative parameter updates. The iterative process can be performed according to... Figure 4 The curves shown are recorded in the following way.

[0154] After a certain bifurcation parameter update, the parent node's radius attribute is 1.20 mm, and the two child node radii are 0.90 mm and 0.80 mm respectively. The sum of the squares of the child node radii is 1.45, which is greater than the sum of the squares of the parent node's radius (1.44), violating the diameter conservation constraint. The electronic device corrects the child node radius attributes by scaling and regenerates the action sequence and predicted structure. Finally, the electronic device calculates the displacement increment, coverage area increment, and volume increment on the calibrated predicted trajectory and combines them to obtain the scalability index. The scalability index output also carries weighting coefficients, residuals, missing rate, break completion records, and iteration records for traceability.

[0155] It should be noted that the embodiments of the present invention have better implementability and are not intended to limit the present invention in any way. Any person skilled in the art may use the above-disclosed technical content to change or modify it into equivalent effective embodiments. However, any modifications or equivalent changes and modifications made to the above embodiments based on the technical essence of the present invention without departing from the content of the technical solution of the present invention shall still fall within the scope of the technical solution of the present invention.

Claims

1. A method for quantitatively evaluating the extensibility of rhizomatous plants, characterized in that, Performed by an electronic device, including the following steps: S1, acquire root and stem imaging data and environmental variable data of the same rhizomatous plant to be evaluated at at least two sampling times; perform three-dimensional reconstruction on the root and stem imaging data at each sampling time to obtain the observation point set. and Rhizome skeleton diagrams were generated based on the observation point set. and ,in, The sampling time number, A set of skeleton nodes. Set up a skeleton edge set; configure node attributes for the skeleton nodes of the rhizome skeleton graph and configure edge attributes for the skeleton edges. The node attributes include at least node coordinates and node radius attributes, and the edge attributes include at least edge length attributes. S2, from the rhizome skeleton diagram Determine the set of extended frontier nodes And construct a digital twin model of the rhizome, defining the model state as Define the environment input as Define the model parameters as a parameter vector. ; S3, Establish the set of extended frontier nodes The set of leading-edge actions, which includes at least: extension actions, bifurcation actions, turning actions, and termination actions; based on the environmental input... With the parameter vector For the extended frontier node set Generate a leading-edge action sequence and apply the leading-edge action sequence to the rhizome skeleton diagram. Obtain the predicted skeleton map ; S4, Construct the consistency objective function The consistency objective function Includes: the set of observation points With the predicted skeleton map Geometric consistency term obtained from the generated prediction point set The rhizome skeleton diagram With the predicted skeleton diagram Topology consistency terms obtained from differences in topology editing actions And the phenotypic consistency term derived from the difference between the observed phenotypic vector and the predicted phenotypic vector. And satisfy: ， in, , These are the weighting coefficients; S5, determine the confidence index based on the observation data quality and registration stability, and adjust the weight coefficients according to the confidence index. , Perform adaptive determination; S6, based on the consensus objective function For the parameter vector Perform iterative calibration, and in each iteration, update the parameter vector. Perform constraint projection to satisfy the constraint set The constraint set It includes at least: minimum spacing constraints between skeleton nodes, bifurcation angle constraints, and diameter conservation constraints; and is based on the calibrated parameter vector. Generate predicted state trajectories; S7, Calculate the scalability index based on the predicted state trajectory. And output the evaluation results, wherein the scalability index At least with the extended frontier node set The displacement increment, coverage area increment, and volume increment are related.

2. The method as described in claim 1, characterized in that, In step S2, the extended frontier node set By screening the rhizome skeleton diagram The terminal node with a degree of 1 is obtained, and the root node corresponding to the starting point of the rhizome is removed from the set of the expanding front node. Excluded from the list.

3. The method as described in claim 2, characterized in that, The end node that intersects with the preset imaging space boundary is marked as the boundary front node, and when generating the front movement sequence in step S3, boundary constraints are applied to the boundary front node to limit its movement type or movement parameters.

4. The method as described in claim 1, characterized in that, In step S3, the set of leading edge actions includes: an extension action that updates the position of the leading edge node along the extension direction, a forking action that generates at least one sub-branch of the leading edge node, a turning action that updates the extension direction of the leading edge node, and a termination action that stops the growth of the leading edge node.

5. The method as described in claim 4, characterized in that, The forward action sequence is represented in the form of a triple of "forward node - action type - action parameter", wherein: the action parameters of the extension action include at least the extension step size and extension direction, the action parameters of the forking action include at least the number of sub-branches and the set of sub-branch directions, the action parameters of the turning action include at least the direction rotation angle, and the action parameters of the termination action include at least the termination marker and the termination time index.

6. The method as described in claim 5, characterized in that, The forking action is defined by a forking rule base, which includes: the conditions for allowing forking at the leading edge node, the range of possible forking child nodes, and the range of possible initial directions for the child branches; and at least a portion of the rule parameters in the forking rule base are defined by the parameter vector. Provided.

7. The method as described in claim 1, characterized in that, In step S4, the topology consistency item By the rhizome skeleton diagram With the predicted skeleton diagram These are converted into topology editing action sequences and aligned to obtain the results. The topology editing actions include at least: node insertion, node deletion, simultaneous insertion, simultaneous deletion, and node attribute replacement; The node attribute replacement includes at least the replacement of the node radius attribute and the edge length attribute; and the topology editing action cost is determined based on the action type cost and the attribute difference cost, and the topology editing action sequence corresponding to the minimum total cost is used to determine the topology consistency item. .

8. The method as described in claim 1, characterized in that, In step S6, the diameter conservation constraint stipulates that the sum of the squares of the radii of the multiple sub-branches obtained by bifurcation at the same parent node does not exceed the square of the radius of the parent node; the bifurcation angle constraint stipulates that the angle between the newly added skeleton edge and its parent skeleton edge falls within a preset angle range; the minimum spacing constraint of skeleton nodes stipulates that the distance between any two different skeleton nodes is not less than the preset minimum spacing.

9. The method as described in claim 8, characterized in that, The parameter vector Including continuous parameters and discrete regular parameters, the constraint projection includes at least: performing boundary clipping on continuous parameters, performing nearest-neighbor replacement within the feasible value set on discrete regular parameters, and performing correction on bifurcation action parameters that violate diameter conservation constraints; and regenerating the frontier action sequence after completing the constraint projection to update the predicted skeleton graph. .

10. The method as described in claim 9, characterized in that, The predicted state trajectory is generated using a rolling time-domain method. When root and stem imaging data at a new sampling time is obtained, the parameter vector is adjusted accordingly. Perform backtracking calibration on the values ​​at least two most recent sampling times, and then backtrack calibrate the parameter vector. The leading-edge action sequence is used for subsequent time-domain generation; and the rhizome skeleton diagram and the calibrated parameter vector are used for this purpose. The consensus objective function The iteration record, the scalability index Write to the storage unit to generate a structured evaluation report.

Images