Bayesian hierarchical model-based interval-valued functional data regression prediction method and system

Abstract

The present application relates to big data processing and meteorological environment monitoring technical field, specifically relates to a kind of interval value function type data regression prediction method and system based on bayesian hierarchical model, comprising: collecting interval value function type observation data;By constructing the function type data set containing interval value covariate, provide accurate input data for regression model;Based on bayesian hierarchical regression model, define and configure the error covariance model and prior model of time-varying binary Gaussian process, generate corresponding error covariance parameters and prior parameters;Based on regression model, construct precision matrix, and use the mixed Markov chain sampling framework of Kronecker product structure to decompose precision matrix, carry out the iterative update of parameter, until parameter converges, output posterior sample.The present application realizes the premise without simplifying model assumption, simultaneously realizes to improve calculation efficiency and prediction accuracy, satisfies the core demand of large-scale complex data analysis.

Description

Technical Field

[0001] This invention relates to the fields of big data processing and meteorological environment monitoring technology, specifically to a method and system for interval-valued function-type data regression prediction based on a Bayesian hierarchical model. Background Technology

[0002] With the advancement of digitalization, fields such as high-frequency financial trading and complex environmental monitoring have generated massive amounts of data. To address the challenges of storage and analysis, interval value data has emerged as an efficient form of data aggregation. It summarizes micro-observations through closed intervals, preserving the variability of the data while capturing measurement errors or inherent fluctuations. For example, in meteorology, the interval formed by the daily maximum and minimum temperatures reflects the overall picture of temperature changes better than the daily average temperature; in financial markets, the interval formed by the intraday maximum and minimum prices vividly portrays market volatility and risk. As the sampling frequency increases, the observed objects often exhibit a close interweaving of interval structure and function structure, i.e., interval value function-type data. For example, daily temperature is not only an interval of maximum and minimum values, but its intraday variation also forms a continuously evolving curve. However, regression analysis for this type of data faces many challenges. Existing technical models make overly simplistic assumptions: the traditional center method only uses midpoint information, ignoring the volatility implied by the interval length; while the center-range method considers the endpoints, it often assumes that the center and range are independent, which is inconsistent with reality. Parameterization methods map intervals to scalars through weighted summation, offering high flexibility. However, in high-dimensional basis function expansions, they are highly susceptible to severe multicollinearity problems, leading to unstable least-squares estimations. Computational efficiency is a significant bottleneck: Bayesian inference faces enormous computational pressure in large-scale functional data scenarios. For data containing… Each observation object and For a dataset with n time sampling points, traditional Bayesian regression algorithms require a dataset with a dimension of n. The global covariance matrix is ​​inverted. This operation has a time complexity of up to [missing information]. As the number of observation stations or the sampling frequency increases, the computational load explodes cubically. For example, in meteorological monitoring networks, and These numbers can reach thousands, making it almost impossible for traditional algorithms to complete the computation within a reasonable time under current hardware conditions, thus failing to meet the timeliness requirements of real-time monitoring and risk warning.

[0003] Furthermore, existing technologies often overlook heteroscedasticity: current interval regression models typically assume homoscedasticity, ignoring the characteristics of data volatility over time. For example, the uncertainty in temperature forecasts often differs significantly before and after sunrise from the afternoon's high temperatures. Ignoring this leads to a significant deviation of forecast interval coverage from the nominal level, underestimating forecast risk. Summary of the Invention

[0004] The purpose of this invention is to provide a method and system for interval-valued function-based data regression prediction based on a Bayesian hierarchical model to solve the above-mentioned technical problems. This method aims to simultaneously improve computational efficiency and prediction accuracy without simplifying model assumptions, thus meeting the core requirements of large-scale complex data analysis.

[0005] The objective of this invention can be achieved through the following technical solutions: A regression prediction method for interval-valued functional data based on a Bayesian hierarchical model is disclosed. The method includes: S1, collecting interval-valued functional observation data, using the interval-valued function curve as the model response variable, constructing a functional dataset containing interval-valued covariates, and dividing the dataset into independent observation units; S2, constructing a Bayesian hierarchical interval-valued function regression model, defining the weight parameters of the parameterization method and the regression coefficients used to capture dynamic dependencies, configuring the error covariance model and prior model of a time-varying bivariate Gaussian process, and generating the corresponding error covariance parameters and prior parameters; S3, constructing an accuracy matrix based on the regression model, and decomposing the accuracy matrix using a hybrid Markov chain sampling framework with a Kronecker product structure; initializing the sampling parameters of the sampling framework and constructing a basis function matrix; under the constraints of the prior parameters, iteratively updating the weight parameters, regression coefficients, and error covariance parameters, reducing the complexity of high-dimensional matrix operations through matrix decomposition, until the weight parameters, regression coefficients, and error covariance parameters converge, and outputting posterior samples including the weight parameters, regression coefficients, and error covariance parameters; S4: Calculate the prediction result of the response variable based on the posterior sample.

[0006] Furthermore, step S1 includes: acquiring the upper and lower bounds of the monitoring object's interval, and using the interval value function curve as the model's response variable. ,in The index of the independent observation unit, and and These represent the lower and upper bounds of the response variable, respectively; the response variable is normalized and expanded using B-spline basis functions to construct a system containing... A functional dataset of observed objects.

[0007] Furthermore, step S2 includes: defining the Bayesian hierarchical interval-valued function regression model expression as: in, This is a function of the lower or upper bound of the response variable. As the reference intercept function, For the first Each covariate on the boundary The regression coefficients, Define the model error term; define the weight parameters of the parameterization method. The formula is: in, and The first Lower and upper bounds of the interval covariates, The error covariance parameter is set to follow a time-varying bivariate Gaussian process, based on logarithmic B-splines. The error covariance model is established using the following formula: in, The correlation coefficient between the upper and lower bounds. The time-varying variance function and error covariance model are used for modeling; the prior model is established based on Bayesian P-splines, and a smoothing parameter is introduced. and auxiliary variables The regression coefficients are subject to structured ridge regularization constraints, as shown in the formula: .

[0008] Furthermore, step S3 is specifically divided into: S31, using a Markov chain sampling framework accelerated by Kronecker product, initializing the sampling parameters of the sampling framework, and constructing the basis function matrix. S32. Transform high-dimensional matrix operations into low-dimensional matrix operations using Kronecker integral solutions, and alternately update the weight parameters and regression coefficients based on random sampling: update the weight parameters λ in the marginal posterior distribution using a collapsed Gibbs strategy; based on the updated weight parameters λ, accelerate the construction of the accuracy matrix using Kronecker product and update the regression coefficients β; update the heteroscedasticity parameters using a Langevin-type sampling method; S33. Determine whether the preset number of sampling steps has been reached based on the current iteration step number; if not, return to S32 and re-execute; if so, perform convergence diagnosis on the Markov chains corresponding to the weight parameters, regression coefficients, and error covariance parameters, and retain valid posterior samples.

[0009] Furthermore, in step S32, the regression coefficients are updated by solving a system of sparse triangular equations.

[0010] Furthermore, step S31 includes: constructing an object-level design matrix. With basis function matrix The precision matrix is ​​obtained by utilizing the Kronecker product property. Decomposition, the formula is: Synchronous initialization of regularization penalty matrix The inverse of the error covariance model The inverse of the error covariance model is a block diagonal sparse matrix.

[0011] Furthermore, step S32 includes: S321, Read the... In the parameter state of the round of iterations, the collapsed Gibbs sampling step is performed: the high-dimensional regression coefficients are obtained through analytical integration. Marginalization from the joint posterior, constructing information about the weight parameters. Marginal log-posterior probability density; S322, in the In each round of iteration, for each weight parameter Perform Metropolis-Hastings update: Propose candidate values ​​through random walks in the Logit transformation space and calculate the acceptance rate. The acceptance rate The calculation includes the marginal likelihood function ratio, the prior probability ratio, and the Jacobian determinant adjustment term introduced by the Logit transformation; wherein, the marginal likelihood function is calculated quickly using the Kronecker product structure; S323, Based on the updated weight parameters Construct the updated object-level design matrix Calculate the conditional posterior accuracy matrix of the regression coefficient β. And introduce a small perturbation term on the diagonal. The formula is: in, This is the error accuracy matrix for a single observation unit. The identity matrix; the error covariance parameter includes the basis coefficients. and error parameters The base coefficient The error parameter Used to describe the correlation between error terms; S324. In this iteration, the B-spline basis matrix is ​​used. Calculating the low-dimensional sparse core matrix of a banded sparse structure Assembly accuracy matrix and perform Cholesky decomposition. ; S325. Sampling regression coefficients by solving the trigonometric equation system The updated formula is: Solving the lower trigonometric equations yields the intermediate vector. And solve the system of upper triangular equations: Output the regression coefficients; where, It is a standard normal random vector; S326. Update the error covariance parameter: Update the error parameter using Fisher-Z transform. The MALA algorithm is used to update the basis coefficients of the heteroscedasticity function using log-posterior gradient information. The formula is: Where ϵ represents the step size and p represents the posterior gradient information; S327. Save the parameter set obtained in the current round. As a posterior sample, it proceeds to the next iteration.

[0012] Furthermore, step S33 includes: Set the number of iterations during the pre-burn-in period and the total number of iterations; during the pre-burn-in period, dynamically adjust the step size in step S326 based on the acceptance rate. To achieve the target acceptance rate; after the pre-burn period, S321 to S327 are executed with a fixed step size, and the samples are stored sparsly at preset intervals to reduce autocorrelation.

[0013] Furthermore, step S4 includes: Calculate the mean of the posterior samples retained after convergence in step S3, and use it as the optimal model parameter estimate. The prediction accuracy of the model was verified using the integral mean square prediction error index; based on the validated optimal weight parameters... and regression coefficients Substitute the interval-valued covariate data for the period to be predicted into the regression model constructed in step S2, and calculate and output the interval upper and lower bound prediction curves of the target response variable. .

[0014] A Bayesian hierarchical model-based interval-valued function-based data regression prediction system includes: a memory and a processor; the memory stores a computer program; the processor executes the computer program to implement the above-described method; the processor is configured to execute the following modules: A data acquisition module for real-time acquisition of interval-valued functional observation data generated by multi-source sensors; a memory management module for constructing and storing the accuracy matrix of regression models; a hybrid sampling execution module for executing a hybrid Markov chain sampling framework with Kronecker product structure; a heteroscedasticity update module for iteratively updating the weight parameters, regression coefficients, and error covariance parameters; and a prediction output module for calculating the prediction results of the response variable based on the posterior samples.

[0015] Compared with the prior art, the present invention has the following beneficial effects: 1. This invention constructs a Bayesian hierarchical interval-valued function regression model and configures an error covariance model for a time-varying bivariate Gaussian process. Based on dynamic dependencies, it effectively estimates the weight parameters and regression coefficients in the regression model. Furthermore, by generating prior parameters and error covariance parameters, it improves the predictive ability of the model.

[0016] 2. This invention achieves the dimensionality reduction of high-dimensional matrix operations by constructing a precision matrix and employing a hybrid Markov chain sampling framework with a Kronecker product structure. Furthermore, it reduces computational complexity by decomposing the precision matrix, thereby accelerating the parameter update and iteration process. Through iterative updates of weight parameters, regression coefficients, and error covariance parameters within the sampling framework, efficient parameter estimation is achieved. Convergence criterion ensures the reliability and stability of the final sampling results, providing highly accurate posterior samples. Thus, it simultaneously improves computational efficiency and prediction accuracy without simplifying model assumptions, meeting the core requirements of large-scale complex data analysis. Attached Figure Description

[0017] Figure 1 This is a flowchart of an interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to the present invention; Figure 2 The baseline trend function estimated from Irish meteorological data in an embodiment of the present invention is shown below. and dynamic heteroscedasticity function Line graph; Figure 3 The dynamic regression coefficient function of each covariate estimated using Irish meteorological data in an embodiment of the present invention is shown below. Line graph. Detailed Implementation

[0018] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments. These embodiments are based on the technical solution of the present invention and provide detailed implementation methods and specific operating procedures. However, the scope of protection of the present invention is not limited to the following embodiments.

[0019] like Figures 1 to 3 The method described is a regression prediction method for interval-valued function data based on a Bayesian hierarchical model. The method includes: S1, collecting interval-valued function observation data, using the interval-valued function curve as the model response variable, constructing a functional dataset containing interval-valued covariates, and dividing the dataset into independent observation units; S2, constructing a Bayesian hierarchical interval-valued function regression model, defining the weight parameters of the parameterization method and the regression coefficients used to capture dynamic dependencies, configuring the error covariance model and prior model of the time-varying bivariate Gaussian process, and generating the corresponding error covariance parameters and prior parameters; S3, constructing an accuracy matrix based on the regression model, and decomposing the accuracy matrix using a hybrid Markov chain sampling framework with a Kronecker product structure; initializing the sampling parameters of the sampling framework and constructing the basis function matrix; under the constraint of prior parameters, iteratively updating the weight parameters, regression coefficients, and error covariance parameters, reducing the computational complexity of high-dimensional matrices through matrix decomposition, until the weight parameters, regression coefficients, and error covariance parameters converge, and outputting posterior samples including the weight parameters, regression coefficients, and error covariance parameters; S4: Calculate the prediction results of the response variable based on the posterior sample.

[0020] Specifically, step S1 includes: collecting interval-valued function type observation data and constructing a dataset. This embodiment first obtains hourly observation records from 24 automatic weather stations of the Met Éireann from January 1, 2018 to December 31, 2022. For the very few missing records caused by sensor malfunctions in the original temperature data, linear interpolation is used to fill in the gaps. In the response variable construction stage, the 24-hour time domain is normalized, and m=51 equidistant grid points are set as the observation benchmark; at each observation time point t, the highest and lowest temperatures observed at all 24 stations across the country are extracted as extreme value boundaries, thereby generating an interval-valued function curve response variable reflecting the daily temperature evolution process across the country. In the covariate construction phase, relative humidity, sea-level pressure, and wind speed observed during the same period were selected. By extracting the global minimum and maximum values ​​of each physical quantity over a 24-hour period each day, interval-valued scalar covariates were formed. Finally, the processed data was divided into independent observation units according to monitoring stations and dates, resulting in a total sample size of [missing data]. The complete dataset was used. To eliminate the influence of differences in the dimensions of physical quantities such as relative humidity, sea level pressure, and wind speed on the calculation, Z-score standardization was performed on all covariates to ensure that the upper and lower bounds of all intervals were within similar orders of magnitude, thereby significantly improving the convergence speed of the Markov chain in the subsequent hybrid MCMC sampling algorithm.

[0021] This embodiment leverages the algebraic properties of the Kronecker product to transform the inversion of high-dimensional dense matrices in Bayesian inference into low-dimensional sparse matrix operations, significantly reducing computational complexity and enabling efficient inference for large sample data. By introducing a hierarchical Bayesian P-spline prior as a structured ridge regularizer, it smooths function curves and suppresses multicollinearity caused by high-dimensional parameter spaces, improving the stability of parameter estimation. Through a time-varying bivariate Gaussian process, it simultaneously captures the correlation between the upper and lower bounds of the interval and intraday dynamic heteroscedasticity, significantly improving the coverage and accuracy of the prediction interval compared to traditional methods that assume homoscedasticity. Utilizing a parameterized method, it automatically learns weights without requiring manual setting of interval feature points and can automatically identify key influencing factors.

[0022] The principle behind the computational efficiency optimization of the algorithm in this invention is analyzed as follows: In the standard Bayesian inference process, the main computational bottleneck lies in the regression coefficients. The posterior distribution is updated. The update process requires constructing the conditional posterior accuracy matrix. For large-scale functional data, the total number of observation points The values ​​are often very large, making it difficult for traditional algorithms to handle their computational load. In traditional regression algorithms, if the accuracy matrix is ​​directly processed... Its dimensions reached Furthermore, the process of constructing the accuracy matrix involves the error covariance parameter matrix. The inverse operation. Because... The dimension is Without utilizing any matrix structure properties, the time complexity of performing the inversion operation will reach [value missing]. This can cause computational crashes when processing high-frequency or large-sample data.

[0023] This embodiment accelerates the process using the following tensor product algebraic identity: Combined with the construction of object-level design matrices With basis function matrix The core terms of the precision matrix can be derived as follows: As can be seen from the above mathematical derivation, this embodiment will process the original process that needed to be performed. Operations on gigantic dense matrices of order 1 can be decomposed into... object-level matrix of order and basis function matrix of order Independent operations. This transformation enables the construction of the precision matrix. The computational complexity was reduced from the original cubic order to Due to the number of covariates With the number of basis functions In practical applications, these are usually small, fixed constants, and the algorithm's time complexity is adjusted based on the sample size. This represents a leap from the cubic to the linear order. This algebraic-level dimensionality reduction transformation forms the theoretical core of this embodiment for processing ultra-large-scale interval-valued function datasets, thereby improving computational efficiency.

[0024] Specifically, step S2 includes: Defining the Bayesian hierarchical interval-valued function regression model expression as: in, This is a function of the lower or upper bound of the response variable. As the reference intercept function, For the first Each covariate on the boundary The regression coefficients, This is the model error term; Define the weight parameters of the parameterized method The formula is: in, and The first Lower and upper bounds of the interval covariates, .

[0025] The error covariance parameter is set to follow a time-varying bivariate Gaussian process, based on logarithmic B-splines. The error covariance model is established using the following formula: in, The correlation coefficient between the upper and lower bounds. The time-varying variance function and error covariance model are used for modeling; a prior model is established based on Bayesian P-splines, and a smoothing parameter is introduced. and auxiliary variables The structured ridge regularization constraint is applied to the regression coefficients, and the formula is as follows: .

[0026] In step S2, the settings are as follows: covariates For humidity, For air pressure, Let the wind speed be denoted by the number of cubic B-spline basis functions. Define the error covariance structure, allowing the variance to vary over time, and specify that the correlation coefficient follows a uniform prior. Parameterize the operator. The weight parameters The prior distribution is set as .

[0027] Specifically, step S3 is divided into: S31, using a Markov chain sampling framework accelerated by Kronecker product, initializing the sampling parameters of the sampling framework, and constructing the basis function matrix. ; Establish the initial state for Markov chain sampling and transform the functional data into a finite-dimensional representation by constructing a basis function matrix; S32, transform high-dimensional matrix operations into low-dimensional matrix operations using Kronecker integral solutions, and alternately update the weight parameters and regression coefficients based on random sampling methods: update the weight parameters λ in the marginal posterior distribution using a collapsed Gibbs strategy; based on the updated weight parameters λ, accelerate the construction of the accuracy matrix and update the regression coefficients β using Kronecker product; update the heteroscedasticity parameters using Langevin-type sampling methods; this is used to reduce the complexity of high-dimensional matrix operations while iteratively updating the weight parameters, regression parameters, and error parameters through block sampling, thereby improving parameter estimation efficiency and enhancing model stability; S33, determine whether the preset number of sampling steps has been reached based on the current iteration step number; if not, return to S32 and re-execute; if so, perform convergence diagnosis on the Markov chains corresponding to the weight parameters, regression coefficients, and error covariance parameters, and retain valid posterior samples. This method is used to improve parameter estimation efficiency and enhance model stability by iteratively updating weight parameters, regression parameters, and error parameters through block sampling while reducing the complexity of high-dimensional matrix operations.

[0028] Specifically, in step S32, the regression coefficients are updated by solving a system of sparse triangular equations. This avoids high-dimensional matrix inversion operations and improves computational efficiency.

[0029] Specifically, step S31 includes: constructing an object-level design matrix. With basis function matrix Using the Kronecker product property to refine the precision matrix Decomposition, the formula is: Synchronous initialization of regularization penalty matrix Inverse of the error covariance model The inverse of the error covariance model is a block diagonal sparse matrix.

[0030] Step S32 includes: S321, Read the... In the parameter state of the round of iterations, the collapsed Gibbs sampling step is performed: the high-dimensional regression coefficients are obtained through analytical integration. Marginalization from the joint posterior involves separating the regression coefficients from the joint posterior distribution, focusing only on the portion related to the weight parameter λ, thus simplifying computation. Through this marginalization process, a posterior distribution related to the weight parameter λ is constructed, which accurately describes the uncertainty of λ. A marginal log-posterior probability density of the weight parameter λ is then constructed. S322, in the In each round of iteration, for each weight parameter Perform Metropolis-Hastings update: Candidate values ​​are proposed through random walks in the Logit transform space, employing nonlinear transformations to adapt to the complex parameter space and increase sampling diversity. Calculate the acceptance rate. Acceptance rate The calculation includes the marginal likelihood ratio, the prior probability ratio, and the Jacobian determinant adjustment term introduced by the Logit transform. The acceptance rate determines the validity of the sampling results based on the marginal likelihood ratio, the prior probability ratio, and the Jacobian adjustment term in the Logit transform, ensuring that the sampling process meets the requirements of Bayesian inference. Specifically, the marginal likelihood function is calculated quickly using a Kronecker product structure. S323, Based on the updated weight parameters Construct the updated object-level design matrix It is used to integrate the relationship between weight parameters and other variables; to calculate the conditional posterior accuracy matrix of the regression coefficient β. And introduce a small perturbation term on the diagonal. Small perturbation term To avoid numerical instability caused by excessively small matrix values ​​or matrices approaching singularity, the formula is: in, This is the error accuracy matrix for a single observation unit. It is the identity matrix; the error covariance parameter includes the basis coefficients. and error parameters base coefficient Error parameters Used to describe the correlation between error terms; S324. In this iteration, the B-spline basis matrix is ​​used. Calculating the low-dimensional sparse core matrix of a banded sparse structure Assembly accuracy matrix B-spline basis matrices have a banded sparse structure, which is used in computation to reduce computational cost, avoid unnecessary calculations, and perform Cholesky decomposition. ; S325. Sampling regression coefficients by solving the trigonometric equation system The regression coefficients are obtained by solving the lower and upper triangular equations; this avoids matrix inversion. The updated formula is: A standard normally distributed random vector z is used as a perturbation to ensure that each update is random; the intermediate vector is obtained by solving the lower trigonometric equations. And solve the system of upper triangular equations: Output the regression coefficients; where, It is a standard normal random vector; S326. Update the error covariance parameter: Update the error parameter using Fisher-Z transform. The MALA algorithm is used to update the basis coefficients of the heteroscedasticity function using log-posterior gradient information. The formula is: Where ϵ represents the step size and p represents the posterior gradient information; S327. Save the parameter set obtained in the current round. As a posterior sample, it proceeds to the next iteration.

[0031] Step S33 includes: setting the number of iterations during the pre-burn-in period and the total number of iterations; and dynamically adjusting the step size in step S326 based on the acceptance rate during the pre-burn-in period. To achieve the target acceptance rate, after the pre-burning period, steps S321 to S327 are executed with a fixed step size, and samples are stored sparsly at preset intervals to reduce autocorrelation. During the sampling process, posterior samples are screened and stored at preset intervals to reduce autocorrelation between samples.

[0032] Building object-level design matrices in memory With basis function matrix . Use Kronecker integral to solve the formula This avoids explicitly storing a huge design matrix. During iterative sampling, a collapsing Gibbs strategy is used to sample and update weights. Wandering in the Logit space and utilizing Calculate the acceptance rate; and use the decomposition formula. Update regression coefficients It only requires inverting a matrix with a dimension of approximately 120×120; the MALA algorithm utilizes gradient information. Updated base coefficients Step length The algorithm was dynamically adjusted to maintain an acceptance rate of approximately 0.574. The algorithm ran for 50,000 iterations, discarding the first 25,000 iterations as a pre-burn period, and then retaining one sample every 5 iterations thereafter, ultimately obtaining 5,000 posterior samples.

[0033] Specifically, step S4 includes: Calculate the mean of the posterior samples retained after convergence in step S3, and use it as the optimal model parameter estimate. The prediction accuracy of the model was verified using the integral mean square prediction error index; based on the validated optimal weight parameters... and regression coefficients Substitute the interval-valued covariate data for the period to be predicted into the regression model constructed in step S2, and calculate and output the interval upper and lower bound prediction curves of the target response variable. .

[0034] To verify the effectiveness of the method in this embodiment, experiments with different sample sizes were designed. and different noise levels Simulation experiments were conducted. Integral mean square prediction error (IMSPE) was used as the prediction accuracy index, and integral mean square error (IMSE) was used as the parameter estimation accuracy index. The comparison results between the method of this embodiment (BHIFR) and the traditional frequency response methods (F-CM, F-CRM, F-MinMax) are shown in Tables 1 and 2 below.

[0035] Table 1. Comparison of out-of-sample prediction performance under different homoscedasticity assumptions

[0036] Table 2 Comparison of parameter estimation errors for benchmark intercept and coefficient of determination

[0037] As shown in Tables 1 and 2, the prediction error and parameter estimation error of the method in this embodiment are significantly lower than those of traditional frequentist methods in all simulation scenarios. This is especially true in challenging high-noise (C=10) and small-sample (n=50) scenarios, where this embodiment demonstrates a significant advantage. This indicates that, compared to traditional frequentist methods that are susceptible to noise overfitting, this embodiment effectively suppresses noise overfitting through the regularization effect of Bayesian hierarchical priors, while maintaining the model's flexibility and prediction robustness. Particularly for zero-effect coefficients... This embodiment demonstrates extremely low estimation error, confirming the model's powerful capabilities in signal recognition and noise reduction, effectively suppressing spurious fluctuations caused by irrelevant variables. The posterior parameter estimation results based on Irish meteorological data reveal the physical mechanism. Humidity weighting. (Indicating that it is mainly affected by the lowest humidity), air pressure Wind speed This indicates that the correlation is mainly influenced by the maximum value, revealing a significant asymmetric dependency mechanism. (Correlation coefficient) This confirms the strong correlation between the upper and lower bounds. For example... Figure 3 As shown, the estimated standard deviation function The model exhibits a significant "bimodal" structure, with the largest fluctuations occurring at sunrise and afternoon, indicating that the model can accurately capture dynamic heteroscedasticity characteristics.

[0038] To verify the adaptability of this embodiment to large-scale data, the physical runtime (in seconds) required for the algorithm of this embodiment and the unoptimized standard algorithm to complete the same number of iterations were recorded under different sample sizes (N), and the speedup ratio was calculated. The results are shown in Table 3 below:

[0039] As shown in Table 3, the computation time of the standard algorithm increases cubically with the increase of the sample size N, which becomes a bottleneck in high-dimensional scenarios. In contrast, the algorithm in this embodiment successfully reduces computational complexity by fully utilizing the sparse Kronecker product structure of the design matrix, demonstrating excellent linear scalability. Specifically, in the case of a large-scale N=1000, this embodiment can complete the computation in only 1.63 seconds, achieving a speedup of over 800 times compared to the standard algorithm, strongly verifying the efficiency and practicality of the algorithm in processing large-scale functional data.

[0040] To further verify the stability of the method in this embodiment, a detailed sensitivity analysis was also conducted on the key hyperparameters affecting model performance. The analysis results are shown in Table 4. Regarding the basis function dimension... The selection of [the option / method] was tested in this embodiment. Different settings. Experimental results show that when At this point, the model achieves an optimal balance between its ability to capture complex curve features and its computational cost. When the dimension is further increased to 20, although the parameter space dimension increases, thanks to the hierarchical Bayesian P-spline prior introduced in this embodiment as a structured ridge regularizer, the model can automatically suppress unnecessary fluctuations. The improvement in prediction error IMSPE is less than 0.5%, and no overfitting occurs. This proves that the choice of function expansion dimension in this embodiment has extremely strong robustness. Secondly, regarding the selection of the hierarchical prior Half-Cauchy scale, this embodiment has a more robust prior structure. The scale parameters were tested in the middle. exist Variations within the range. Experiments show that even when the prior distribution changes from a strongly contracted state to a highly diffuse state, the weight parameters remain unchanged. The posterior mean fluctuations are extremely small, with deviations all less than 0.001. This result strongly demonstrates that the hierarchical prior structure adopted in this embodiment has excellent adaptability, effectively extracting smoothing information from the observed data without relying on manually set specific hyperparameter values, thus ensuring the objectivity of the inference results. Finally, regarding the correction measures for computational stability, this embodiment constructs the posterior accuracy matrix in step S323. At that time, a small perturbation term was introduced on its diagonal. Comparative tests revealed that the perturbation term exist to When the value is between these ranges, the impact on the deviation of the final prediction result is almost negligible. However, the introduction of this term greatly improves the condition number of the high-dimensional precision matrix, ensures the numerical stability of Cholesky decomposition under various extreme noise interferences, avoids algorithm collapse caused by matrix singularity, and provides a solid numerical guarantee for processing large-scale real-time data.

[0041] Table 4. Sensitivity Analysis Results of Key Hyperparameters

[0042] A Bayesian hierarchical model-based interval-valued function data regression prediction system includes: a memory and a processor; the memory stores a computer program; the processor executes the computer program to implement the above method; the processor is configured to execute the following modules: A data acquisition module for real-time acquisition of interval-valued functional observation data generated by multi-source sensors; a memory management module for constructing and storing the accuracy matrix of regression models; a hybrid sampling execution module for executing a hybrid Markov chain sampling framework with Kronecker product structure; a heteroscedasticity update module for iteratively updating weight parameters, regression coefficients, and error covariance parameters; and a prediction output module for calculating the prediction results of the response variable based on posterior samples.

[0043] It should be noted that the terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the exemplary embodiments according to this application. As used herein, the singular form is intended to include the plural form as well, unless the context clearly indicates otherwise. Furthermore, it should be understood that when the terms "comprising" and / or "including" are used in this specification, they indicate the presence of features, steps, operations, devices, components, and / or combinations thereof.

[0044] It should be noted that the terms "first," "second," etc., used in the specification, claims, and accompanying drawings of this application are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be interchanged where appropriate so that the embodiments of this application described herein can be implemented in sequences other than those illustrated or described herein.

[0045] The present invention has been further described above with reference to specific embodiments. However, it should be understood that the specific description herein should not be construed as limiting the nature and scope of the present invention. Various modifications made to the above embodiments by those skilled in the art after reading this specification are all within the scope of protection of the present invention.

Claims

1. A method for interval-valued function-based data regression prediction based on a Bayesian hierarchical model, characterized in that, The method includes: S1. Collecting interval-valued function observation data, using the interval-valued function curve as the model response variable, constructing a function-type dataset containing interval-valued covariates, and dividing the dataset into independent observation units; S2. Constructing a Bayesian hierarchical interval-valued function regression model, defining the weight parameters of the parameterization method and the regression coefficients used to capture dynamic dependencies, configuring the error covariance model and prior model of the time-varying binary Gaussian process, and generating the corresponding error covariance parameters and prior parameters; S3. Constructing an accuracy matrix based on the regression model, and decomposing the accuracy matrix using a hybrid Markov chain sampling framework with a Kronecker product structure; initializing the sampling parameters of the sampling framework and constructing a basis function matrix; under the constraints of the prior parameters, iteratively updating the weight parameters, regression coefficients, and error covariance parameters, reducing the complexity of high-dimensional matrix operations through matrix decomposition, until the weight parameters, regression coefficients, and error covariance parameters converge, and outputting posterior samples including weight parameters, regression coefficients, and error covariance parameters; S4: Calculate the prediction result of the response variable based on the posterior sample.

2. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 1, characterized in that, Step S1 includes: acquiring the upper and lower bounds of the monitoring object's interval, and using the interval value function curve as the model's response variable. ,in The index of the independent observation unit, and and These represent the lower and upper bounds of the response variable, respectively; the response variable is normalized and expanded using B-spline basis functions to construct a system containing... A functional dataset of observed objects.

3. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 1, characterized in that, Step S2 includes: Defining the Bayesian hierarchical interval-valued function regression model expression as: in, This is a function of the lower or upper bound of the response variable. As the reference intercept function, For the first Each covariate on the boundary The regression coefficients, Define the model error term; define the weight parameters of the parameterization method. The formula is: in, and The first Lower and upper bounds of the interval covariates, The error covariance parameter is set to follow a time-varying bivariate Gaussian process, based on logarithmic B-splines. The error covariance model is established using the following formula: in, The correlation coefficient between the upper and lower bounds. The time-varying variance function and error covariance model are used for modeling; the prior model is established based on Bayesian P-splines, and a smoothing parameter is introduced. and auxiliary variables The regression coefficients are subject to structured ridge regularization constraints, as shown in the formula: 。 4. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 3, characterized in that, Step S3 is specifically divided into: S31, a Markov chain sampling framework based on Kronecker product acceleration, initializing the sampling parameters of the sampling framework, and constructing the basis function matrix. S32. Transform high-dimensional matrix operations into low-dimensional matrix operations using Kronecker integral solutions, and alternately update the weight parameters and regression coefficients based on random sampling: update the weight parameters λ in the marginal posterior distribution using a collapsed Gibbs strategy; based on the updated weight parameters λ, accelerate the construction of the accuracy matrix using Kronecker product and update the regression coefficients β; update the heteroscedasticity parameters using a Langevin-type sampling method; S33. Determine whether the preset number of sampling steps has been reached based on the current iteration step number; if not, return to S32 and re-execute; if so, perform convergence diagnosis on the Markov chains corresponding to the weight parameters, regression coefficients, and error covariance parameters, and retain valid posterior samples.

5. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 4, characterized in that, In step S32, the regression coefficients are updated by solving a system of sparse triangular equations.

6. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 4, characterized in that, Step S31 includes: constructing an object-level design matrix. With basis function matrix The precision matrix is ​​obtained by utilizing the Kronecker product property. Decomposition, the formula is: Synchronous initialization of regularization penalty matrix The inverse of the error covariance model The inverse of the error covariance model is a block diagonal sparse matrix.

7. The interval-valued function data regression prediction method based on a Bayesian hierarchical model according to claim 6, characterized in that, Step S32 includes: S321, Read the... In the parameter state of the round of iterations, the collapsed Gibbs sampling step is performed: the high-dimensional regression coefficients are obtained through analytical integration. Marginalization from the joint posterior, constructing information about the weight parameters. Marginal log-posterior probability density; S322, in the In each round of iteration, for each weight parameter Perform Metropolis-Hastings update: Propose candidate values ​​through random walks in the Logit transformation space and calculate the acceptance rate. The acceptance rate The calculation includes the marginal likelihood function ratio, the prior probability ratio, and the Jacobian determinant adjustment term introduced by the Logit transformation; wherein, the marginal likelihood function is calculated quickly using the Kronecker product structure; S323, Based on the updated weight parameters Construct the updated object-level design matrix Calculate the conditional posterior accuracy matrix of the regression coefficient β. And introduce a small perturbation term on the diagonal. The formula is: in, This represents the error accuracy matrix for a single observation unit. The identity matrix; the error covariance parameter includes the basis coefficients. and error parameters The base coefficient The error parameter Used to describe the correlation between error terms; S324. In this iteration, the B-spline basis matrix is ​​used. Calculating the low-dimensional sparse core matrix of a banded sparse structure Assembly accuracy matrix and perform Cholesky decomposition. ; S325. Sampling regression coefficients by solving the trigonometric equation system The updated formula is: Solving the lower trigonometric equations yields the intermediate vector. And solve the system of upper triangular equations: Output the regression coefficients; where, It is a standard normal random vector; S326. Update the error covariance parameter: Update the error parameter using Fisher-Z transform. The MALA algorithm is used to update the basis coefficients of the heteroscedasticity function using log-posterior gradient information. The formula is: Where ϵ represents the step size and p represents the posterior gradient information; S327. Save the parameter set obtained in the current round. As a posterior sample, it proceeds to the next iteration.

8. The interval-valued function data regression prediction method based on a Bayesian hierarchical model according to claim 7, characterized in that, Step S33 includes: Set the number of iterations during the pre-burn-in period and the total number of iterations; during the pre-burn-in period, dynamically adjust the step size in step S326 based on the acceptance rate. To achieve the target acceptance rate; after the pre-burn period, S321 to S327 are executed with a fixed step size, and the samples are stored sparsly at preset intervals to reduce autocorrelation.

9. The interval-valued function-based data regression prediction method based on a Bayesian hierarchical model according to claim 1, characterized in that, Step S4 includes: Calculate the mean of the posterior samples retained after convergence in step S3, and use it as the optimal model parameter estimate. The prediction accuracy of the model was verified using the integral mean square prediction error index; based on the validated optimal weight parameters... and regression coefficients Substitute the interval-valued covariate data for the period to be predicted into the regression model constructed in step S2, and calculate and output the interval upper and lower bound prediction curves of the target response variable. .

10. A data regression prediction system based on a Bayesian hierarchical model with interval-valued functions, characterized in that, include: Memory and processor; The memory is used to store a computer program; the processor is used to execute the computer program to implement the method as described in any one of claims 1-9; the processor is configured to execute the following modules: A data acquisition module for real-time acquisition of interval-valued functional observation data generated by multi-source sensors; a memory management module for constructing and storing the accuracy matrix of regression models; a hybrid sampling execution module for executing a hybrid Markov chain sampling framework with Kronecker product structure; a heteroscedasticity update module for iteratively updating the weight parameters, regression coefficients, and error covariance parameters; and a prediction output module for calculating the prediction results of the response variable based on the posterior samples.

Images