Research on the joint fluctuation laws between locational marginal price and renewables based on complex networks: A case study in Independent System Operator New England

The growing share of renewable energy has an increasingly significant impact on the electricity market. Locational marginal price (LMP), as a sensitive price signal, has a potential but considerable correlation with renewable energy generation. In this context, this paper focuses on the joint fluctuation laws between LMP and renewables consumed to generate electricity. At first, LMP and renewables consumption data from Independent System Operator New England (ISO‐NE) are collected and preprocessed. Then, through highly comparative time‐series analysis (hctsa), it is found that discrete symbolization features have effective description ability. Therefore, a coarse‐grained method is conducted to convert the real matrix data into a symbol vector based on numerical distribution analysis. Next, joint fluctuation complex networks are constructed according to the symbol vector. From a global perspective, the joint fluctuation network follows an exponential distribution and exhibits hierarchical modularity. From a local perspective, the joint fluctuation network follows a power‐law distribution. The research results of this paper provide a new perspective to understand the evolution of the electricity market with large‐scale renewable energy and have reference value for identifying the valuable joint fluctuation patterns and improving the LMP forecast results.


| INTRODUCTION
The development and utilization of renewable energy is an inevitable choice for the sustainable development of human society. According to BP Energy Outlook-2019 edition, renewable energy's share of global power will increase from 8.5% in 2018 to 29% by 2040, with solar and wind power playing an important role. 1 Solar and wind energy, as the typical representatives of renewable energy, have already had a significant impact on the energy market and become a factor that cannot be ignored in price decisions.
Locational marginal price (LMP) is generally considered to be a sensitive and effective signal, which is defined as the cost of supplying an increment of load. Grid-connected renewable energy generation will increase the generation randomness. Distributed small-scale renewable power stations will also affect the physical limits of the transmission system. Those factors will have a direct impact on LMP.

| Impact analysis of renewable energy generation on the electricity market
The traditional idea of impact analysis is to consider renewable energy generation as a new factor to improve the accuracy of electricity market forecast. For instance, Rodat et al introduced meteorological factors into the short-term forecast of a solar thermal power plant to improve the accuracy of the daily forecast of solar thermal energy production. 3 In terms of long-term forecast, Yang et al considered the impact of social and economic factors on renewable energy consumption. 4 The trend in such direct forecast studies is the application of modern intelligent algorithms, such as Support Vector Machine (SVM), Grey Model, 5,6 and ANN. 7 Although the accuracy of forecasts continues to improve, it is difficult to explain the impact of renewables on electricity markets due to the uninterpretability of the methods themselves.
Other scholars have broadened research ideas. They explored the equilibrium relationships between renewable energy consumption, power generation costs, and carbon emission, 8 that is, to directly study the relationship between renewable energy and market (both electricity market and carbon market). Among them, the impact of renewables on electricity prices has always been the focus. Zipp 9 studied the price effect of the rising share of renewable electricity by multivariate regression model to conclude the marketability of variable renewable energy in liberalized electricity markets. Hartner and Permoser 10 used dispatch models to analyze the effect of electricity generation from photovoltaic on electricity prices and revenues for storage plants. The methods of these studies are relatively traditional (eg, DEA 8 and VECM 11 ), but they have inspired the research ideas in this paper. If LMP and renewables are regarded as a unified system, their relationship will become nonlinear, dynamic, and complex, and then new methods will be required to reveal the evolution of the system.

| Time-series analysis based on complex networks
As discussed above, new methods are required for multi-factor time-series analysis. Some scholars have explored this and applied interdisciplinary methods such as statistical mechanics and econophysics methods to the analysis with great success.
1. The visibility graph (VG) method. 12,13 Some early studies proposed the method to transform the dynamic linear regression model between two time series into a complex network. 14 On this basis, Gao et al developed a multiscale limited penetrable horizontal visibility graph (MLPHVG) to explore the dynamic behaviors from multiscale analysis and network analysis. 15 MLPHVG is proved to be effective in distinguishing among random, periodic, and chaotic signals. 16 In other words, it retains the periodicity of time series but weakens the sequential relation. 2. The phase space reconstruction method. 17 The idea of phase space reconstruction is to embed the scalar time series into a high dimensional phase space. 18 The basis of the above two methods is the coarse-grained process, which is the core of symbolic dynamics. Symbolic dynamics has unique advantages in time-series analysis since it can simplify the description of time series. Despite the loss of accuracy, some statistical features are more clearly exposed. However, the coarse-grained process is not applicable to all time series. The rationality of coarse-grained process should be demonstrated before a particular study begins. What is more, the design of the coarse-grained formula for a particular time series also affects the reliability of the results. The above interdisciplinary methods have achieved certain application results in many fields. Gao et al applied complex networks to detect crude oil price fluctuation mechanism under different periodic time series. 19 In the financial market, when the research object is extended to two time series, the above methods can be combined with regression analysis for similar studies. 20 Inspired by these studies, this paper adopts symbolic dynamics and complex network theory to study the joint fluctuation laws between LMP and renewables: (a) this paper designs more reasonable and practical analysis methods for specific time series, which makes up for the deficiencies of previous studies. (b) Four time series are considered, extending previous studies. (c) An example verifies the effectiveness and value of the research results.

| Complex network theory
The dominant theory of this paper is complex network theory. Related research on complex network theory has a long history, but modern research originated in 1998 when Watts and Strogatz put forward the small-world network model. 21 Complex network theory is a useful tool to study large-scale networks and also suitable for studying complex evolutional relationships. The statistical indicators used in this paper are as follows:

| Clustering coefficient
The clustering coefficient is a measure of the local cohesiveness of a node. In case of weighted complex networks, let k i be the degree of node n i , ω ij be the weight of the edge between nodes n i and n j , 〈ω i 〉 be the average weight incident to node n i , and a ij be the presence of a connection between nodes n i and n j , and then, the weighted clustering coefficient C i of node n i is defined as 22 The clustering coefficient C of the network is the average of the clustering coefficients of all nodes.

| Modularity
Modularity is a measure of the difference between a network and a random network in a certain community partition, which can indicate the quality of the partition. For a network, different community partition corresponds to different modularity. The larger the modularity, the more reasonable the partition; the smaller the modularity, the more blurred the partition.
In case of weighted networks, let k i = ∑ j ij be the sum of the weights of the edges attached to node n i , and c i be the community to which node n i is assigned. Write m = 1 2 ∑ ij ij . Define the δ function as.
Then the modularity Q of the network is defined as 23

| Modified HITS algorithm
The modified HITS (M-HITS) algorithm is designed based on the HITS algorithm 24,25 to measure the statistical analysis value of each node. The M-HITS algorithm computes two separate values for each node, the out-CRn oc n and the in-CRn ic n . The out-CRn measures the concentration ratio of the edges pointing from the node. And the in-CRn measures the concentration ratio of the edges pointing to the node. These two values reflect the complexity of the joint fluctuation laws from the perspective of concentration.
The algorithm is built on two basic assumptions about mutual enhancement: (a) nodes with high out-CRn will have high-weight edges centrally pointing to the few nodes with high in-CRn and (b) nodes with high in-CRn will have highweight edges centrally pointing from the few nodes with high out-CRn.
Write the number of nodes as n, the out-CRn vector as OC = (oc 1 , oc 2 , …, oc n ), and the in-CRn vector as IC = (ic 1 , ic 2 , …, ic n ). The algorithm steps are described as follows: Step 1: In the initial state, the oc i and the ic i of node n i are both set to 1, wherein, 1 ≤ i ≤ n.
Step 2: For each node, write the weights of all the edges pointing by node n i as W out i = (w out i,1 ,w out i,2 , ⋯ ,w out i,n ) and multiply it by IC element by element to get the vector OA i = W out i . * IC. Similarly, write the weights of all the edges pointing to node n i as W in (1) multiply it by OC element by element to get the vector IA i = W in i . * OC. Note that if there is no edge between two nodes, set the weight to 0.
Step 3: For each OA i and IA i of node n i , sort them in descending order, respectively, to get the vectors SOA i = (soa i1 ,soa i2 , ⋯ ,soa in ) and SIA i = (sia i1 ,sia i2 , ⋯ ,sia in ).
Step 4: For each SOA i of node n i , divide the former term soa i,j by the later term soa i,j+1 and sum up the quotients to get the new oc i . Similarly, for each SIA i of node n i , divide the former term sia i,j by the later term sia i,j+1 and sum up the quotients to get the new ic i . Wherein, 1 ≤ j ≤ n − 1. Note that in order to avoid the exception that the divisor is equal to zero, all weights are increased by one-unit value 1. In mathematical words, update oc i and ic i according to the following formulas: Step 5: Normalize OC and IC through the following formula: Step 6: Repeat Step 2 to Step 5 until the changes of OC and IC are lower than the acceptable accuracy. OC and IC will eventually converge to stable vectors, that is, the final out-CRn and in-CRn.
The M-HITS algorithm flowchart is shown as Figure 1.

| Symbolic dynamics and hctsa
The basic idea of symbolic dynamics is to convert the numeric data into character data (ie, symbols) for analysis. Coarse-grained process is one of the core method of symbolic dynamics. 26 Symbolic dynamics could establish links between the evolution track and the formal language. 27 And the coarse-grained process allows to discard low-level details while retaining high-level features. 28 This paper adopts symbolic dynamics rather than a regression method for three reasons: (a) four factors are considered, thus the regression model will be too complicated and its effect will be too poor when there is less data. A simplified description is required to expose the joint fluctuation laws. (b) The impact of renewables on LMP is limited, and there is no significant correlation between them, so fine regression analysis is not feasible. (c) Highly comparative time-series analysis (hctsa) 29 proves discrete symbolization to be an effective feature for describing the time series of the electricity market with large-scale renewable energy.
hctsa is a MATLAB package that can extract over 7700 features of time series, such as the maximal or minimal value. The features give new insights into time series and their dynamics. hctsa also provides a suite of analysis tools for selecting useful and interpretable features to describe time series or for other applications. It is noteworthy that there is a strong correlation between the wind series and the total renewables series according to Pearson's correlation analysis (the Pearson Correlation Coefficient is .948, and the correlation is significant at .01 level). In addition, hydro is also taken into account because of its large share of energy. The four time series (LMP, hydro, wind, and solar) are shown in Figure 2

| Data preprocessing
The four time series are organized through linear interpolation method to ensure the consistency of statistical time points. Then, the first-order difference of each time series is calculated to eliminate the impact of long-term trends. At last, the nighttime data are removed since there is no photovoltaic output during the night.

| Coarse-grained process
First, to demonstrate the rationality of coarse-grained process, hctsa is implemented on the above four vectors. The performance of the extracted 6330 features is compared using a linear SVM classifier, and the classification accuracy is used to measure the ability to describe the time series. There are 57 features with 100% classification accuracy, including stationarity and step detection, correlation, entropy, and information theory. Among them, there are 14 discrete symbolization features, accounting for 24.6%. In addition, among the top 200 best features, the discrete symbolization features account for 15%. Based on the above analysis, it can be seen that the symbolization method is necessary and effective for further analysis. Combined with the numerical distribution of time series, the following coarse-grained process is adopted and designed to simplify the joint fluctuation patterns.
The core of the coarse-grained process is to determine the granularities. Therefore, the probability density estimate of each vector is plotted as As shown in Figure 3, the four probability density curves all appear as bell curves with clear corners. The curve shape expresses the distribution feature of variable values-most values are concentrated around a certain value. Since the study only focuses on the joint fluctuation laws at a coarse grain level, the granularities can be determined based on the shape of the probability density distribution curve. The corners indicated by the vertical lines in Figure 3 represent the Based on the granularities, time series can be converted into character vectors. Character set {T, R, E, F, B} is adopted to describe the fluctuation of time series, wherein, T indicates a sharp rise, R indicates a rise, E indicates an equal, F indicates a fall, and B indicates a sharp fall. More characters can describe variables more precisely, but increase the number of joint fluctuation patterns, which will perplex the complex networks established later and disguise the joint fluctuation laws. So, five characters are appropriate.
Through the above processing, the real matrix DM is converted into a character matrix. Each row of the character matrix contains four characters, which can be treated as a symbol and represent a joint fluctuation pattern. That is the character matrix can be treated as a symbol column vector MV of length 42 237. The 200 sets of data on the head and tail of MV are excluded for future testing. The first 14 rows of the vector are shown in Figure 4 as samples.

| Joint fluctuation complex networks
Based on the symbol column vector MV, complex networks can be constructed to study the joint fluctuation laws between LMP and renewables. According to the coarsegrained formulas, there will be 5 3 × 3=375 symbols, namely 375 joint fluctuation pattern in theory. But only 324 patterns appear in the study. The 324 patterns N = {n i | i = 1, 2, …, 324} are set as the nodes of the complex networks. Patterns such as TEBR does not appear, which means the joint fluctuation pattern with sharply rising LMP, constant hydro, sharply falling wind, and rising solar does not exist in statistics. The order of the patterns in MV represents the state changes of the system over time. Thus, it is reasonable to consider every two adjacent patterns in MV to represent an edge from the former to the latter. The weight of the edge is 1. And the edge set is written as E = {e i | i = 1, 2, …, 42 236}. By this way, a directed complex network CN = (N, E) is constructed. Figure 5 shows an overview of the joint fluctuation network. Wherein, the color and size of nodes vary with the degrees. That is, the greater the degree, the larger the size and the redder the color. The edges are colored in accordance with the source nodes.

| Global laws of joint fluctuations-the static geometric features of the network
The static geometric features of the joint fluctuation network fall in between regular networks and random networks. The graph density ρ is 0.404 < 0.5, which means that the network is not a complete graph, and the conversions between nodes are not arbitrary. The average path length ℓ is 2.303, which means that any conversion between nodes takes an average path of 2.3 nodes. The diameter of the network Φ is 7, which is the maximum of the shortest paths in the network and represents the extreme conversion case between nodes: It will pass through seven nodes to convert from node TFRE to TREE. It means that if TFRE appears, TREE will not appear in a short term (at least about 37 minutes). In other words, in the pattern TFRE, the case with a significant wind power changes while a stabile system is unlikely to occur. So, it can be guessed that wind power is a more important factor affecting the dynamic behaviors of the system. This guess can be supported by other long paths.

| Degree distribution-the scale-free property of the network
Degree distribution is a core feature of complex networks. Each node in the joint fluctuation network has the same indegree and out-degree except RTFF and REFF. (They are the first and the last pattern of MV, respectively.) So, it is indifference to use the in-degree, the out-degree, or the degree for analysis. This paper uses the degree for analysis. The degree k i of node n i reflects the occurrence probability of the corresponding joint fluctuation pattern, and the degree distribution reflects the concentration or dispersion trend of joint fluctuations. It can be visually seen from Figure 5 that several nodes, such as RRFR, RRRF, and RRFF, have very high degrees, which is a typical feature of scale-free networks. To be specific, there are 127 nodes with degrees greater than the average 〈k〉 = 130.36, and their cumulative degree accounts for 91.34% of the total. The top 33 nodes with the highest degree have a cumulative degree accounting over 50% of the total. These statistics indicate the concentration trend of the network.
The cumulative degree distribution of the joint fluctuation network is shown as Figure 6. Wherein, the horizontal axis represents the degree with the range of [2,2186]; the vertical axis represents the probability; the points (denoted by "×") represent the probability that the degree is greater than the corresponding x value. Note that the vertical axis is a logarithmic axis. The solid lines in Figure 6 represent the best fit to the nodes based on nonlinear least squares method. The expression of the solid line is y = 0.571e −0.002x , R 2 = 0.968. In some literature, similar cases are fitted to power-law distribution to reveal the scale-free feature of the networks. But the fitting effect is sometimes not good enough. [30][31][32] The fitting result in this paper indicates that the degree distribution of the network is deviated from a pure power-law and conformed to an exponential. In practice, several real networks do show similar features. For example, the degree distribution of the electric power grid of southern California is more consistent with a single-scale exponential distribution. 33 Amaral et al attribute the phenomenon to the growth constraints of the network, such as aging and cost. 34 The reason why the degree distribution of the network displays an exponential will be discussed later when analyzing its evolution.

| Clustering-degree correlation-the hierarchical modularity property of the network
Hierarchical modularity is another common feature of complex networks. It is usually manifested as follows: Nodes with low degrees usually have high clustering coefficients and belong to tightly connected small modules. On the contrary, hub nodes with high degrees have low clustering coefficients, and their functions are just to connect different modules. In other words, in a complex network with hierarchical modularity, many small modules with dense internal relations are loosely related to each other, thus forming a larger scale topology module. One of the important marks of hierarchical modularity is that the clustering-degree correlation displays a powerlaw distribution. It should be noted that neither ER random graph nor BA scale-free network has hierarchical topology.
The clustering-degree correlation of the joint fluctuation network is shown in Figure 7. The clustering-degree correlation and the fit lines are plotted in a log-log coordinate. Wherein, the blue "+" describes the relationship between the degrees and the clustering coefficients of nodes, while the red "×" describes the relationship between the degrees of nodes and the average clustering coefficients of their out clusters. 35 (the out cluster of node n i is the small network consisting of the edges and nodes pointed from node n i , and it describes the next evolution of joint fluctuations in this paper).
As can be seen from Figure 7, the clustering-degree correlation displays a power-law distribution from the perspective of the network itself or the out cluster, which indicates that the network is hierarchical modularity. Therefore, it is necessary to detect communities of the network and study its modularity.

| Community detection-the regularity of joint fluctuations
The fast unfolding algorithm proposed by Lambiotte and other scholars 36,37 is adopted to detect communities of the joint fluctuation network in this paper. The algorithm is a heuristic method based on modularity optimization, and its steps are not discussed in this paper. It is adopted for the following three reasons: (a) easy to operate; (b) linear time complexity; (c) a preset number of communities are not required.
The results of community detection are shown in Table  A1. As can be seen from Table A1, the number of communities is six and the modularity is 0.320. In order to observe each community more intuitively, draw radar maps as Figure 8: let CL, CH, CW, and CS be the axes, and place T, R, E, F, and B (or R, E, F) on the axes in order. So, each node can be expressed as a translucent block, and each community can be expressed as a radar map shown as Figure 8(A) to (F).
As can be seen from Figure 8, in community 0, almost all nodes have a CS column of R, and the CW column is centered on B and F; in community 1, the CS column is still always R, but the CW column is centered on R, T, and E; in community 2, the CH column is always T; in community 3, the CS column is always E; in community 4, the CS column is always F; and in community 5, the CH column is always B. In order to further analyze the effectiveness of community detection, it is necessary to research the time-domain distribution of elements (patterns and communities) before and after community detection. Note that the position of elements in vector MV reflects the time-domain distribution, so the function find of MATLAB can be used to obtain the linear indices of elements in MV, which is recorded as vector Ind. Take community 0, community 1, pattern RRFR, and pattern BBRF as examples. Figure 9 shows their time-domain distribution. Wherein, (A) corresponds to community 0, (B) corresponds to community 1, (C) corresponds to pattern RRFR, which is the most frequent pattern, and (D) corresponds to pattern BBRF, which has the most distorted figure. Figure  9(A) and (B) is the result after community detection, (C) and (D) is the result before community detection, it can be visually seen that (C) and (D) are more distorted than (A) and (B). To measure the differences before and after community detection, it is feasible to fit a linear regression equation to Ind and calculate the root-mean-square error (RMSE). Before community detection, there are 324 patterns in MV. A total of 98 patterns are discarded because they appear <10 times and are unsuitable for linear fitting. Based on the reserved patterns, 226 RMSEs are obtained, shown as Figure 10(A). Similarly, there are six communities in MV after community detection, corresponding to 6 RMSEs showing as Figure 10(B). The dimensions of these 226 + 6 sets of data are consistent, both about [1,42237]; thus, it is significative to directly compare the mean of RMSEs. The mean of RMSEs before community detection M bf is 1565.43, and the mean after community detection M aft is 251.39. M aft is significantly less than M bf , which indicates that community detection makes the regularity of joint fluctuations more clearly exposed.
Community mean patterns tend to evolve within specific clusters, which will be further analyzed in Section 5.2.1.

HITS algorithm
Although community detection has made the regularity of joint fluctuations clearer, we still hope to observe the regularity at a more detailed level. It is obviously that nodes with F I G U R E 1 1 Nodes distribution in a out-in-CRn coordinate system F I G U R E 1 2 Cumulative in-degree distribution of the RRFR-out network F I G U R E 1 3 The conversion probability matrix a smaller in or out cluster will have a more concentrated evolution direction, which is of greater statistical analysis value for forecasting and decision making.
The M-HITS algorithm proposed in this paper can measure the statistical analysis values of nodes. According to the M-HITS algorithm, the out-CRn and in-CRn of each node are calculated and shown in Table A1. As can be seen from Table A1, only a few nodes have high statistical analysis values, and they generally have lower degrees (not absolutely). To visualize the results, take out-CRn as the abscissa and in-CRn as the ordinate, plot the nodes in Figure 11. As can be seen from Figure 11, the nodes roughly distributed near the straight line l, which indicates that nodes with high in-CRn generally also have high out-CRn. The lines p and q divide the coordinate system into four regions, namely (a), (b), (c), and (d). Nodes in (c) are considered to have lower statistical analysis values due to their low out-CRn and in-CRn. Nodes in (a) have high in-CRn, nodes in (d) have high out-CRn, and nodes in (a) have both high in-CRn and high out-CRn. These nodes, such as node RRFR, are considered to have higher statistical analysis values and will be further analyzed in the next section.

| Local laws of joint fluctuations-the evolution process of the network
The previous section analyzes the system during a certain period (March 2016-March 2017) from a static perspective. In fact, the system is dynamic, and the network will continue to change with the time window moving.
It was indicated earlier that the generative process of the joint fluctuation network is the cause of the exponential degree distribution. Therefore, it is necessary to study the evolution process of the network. Different from general research, the evolution of the joint fluctuation network is not the addition of new nodes or the removal of old nodes, but the formation of new edges between the last node in MV and other nodes. Therefore, the study on evolution is actually the analysis of the out cluster of a particular node.

| Analysis of the out cluster of a single node
At first, still take node RRFR as an example to explore the laws of its next-step evolution. Record the out cluster of RRFR as the RRFR-out network. The cumulative in-degree distribution of the RRFR-out network is shown in Figure 12. The straight lines in the log-log coordinate system show that the in-degree distribution of the network follows a power-law distribution.
The sample statistics of the RRFR-out network is shown in Table 1. There are 113 nodes in the network. Among them, 25 nodes belong to community 0, the same as node RRFR, accounting for 22.124% of the total, while their cumulative indegree accounts for 69.259%. In addition, except node RRRR, the top eight nodes with the highest in-degree in the network all belong to community 0. All these indicate that the nextstep evolution of node RRFR tends to convert within community 0 in general. To be specific, the next-step evolution of RRFR tends to maintain the XRFR pattern (X represents a character in {T, R, E, F, B}): It has a 32.388% chance to maintain node RRFR, 12.626% to convert to node ERFR and 6.770% to covert to node FRFR. Besides, it has a high probability to convert to node RRRR, RFFR, REFR, EFFR, and RRBR. Similar conclusions can be drawn from the analysis of nodes such as RRRF and RRFF. Based on these calculations, a conversion probability matrix as shown in Figure 13 can be obtained. The element in the matrix represents the transition probability from the source node to the target node and are mapped to the color gamut shown in Figure 13. In the matrix, the sum of the columns is equal to 1. The diagonal elements usually have high values and appear light blue, indicating a certain tendency to remain constant during the evolution.
The above analysis shows that, in the evolution process, new edges do not follow the preferential attachment rules of scale-free networks, but follow the power-law distribution to establish connections between nodes.

| Analysis of the out cluster of a class of nodes
Further, take pattern XRFR as an example to explore the laws of its next-step evolution. The out cluster of pattern XRFR is the network consisting of the edges and nodes pointed from nodes TRFR, RRFR, ERFR, FRFR, and BRFR and is recorded as the XRFR-out network. The overview of the XRFR network is shown as Figure 14, and its cumulative in-degree distribution is shown as Figure 15. The network includes 178 nodes, wherein 26 nodes have an in-degree greater than the average of 16.062, accounting for 14.607% of the total, while their cumulative in-degree accounts for 78.734%. The network tends to evolve among pattern RRFR, ERFR, and FRFR. That is, it tends to keep the pattern of Hydro, Wind and Solar as RFR while fluctuating the LMP in evolution. If these three nodes are excluded, it can be seen that the next-step evolution of XRFR tends to pattern XRRR or XFFR. This may be related to the randomness of the Wind and the compensation and regulation functions of the Hydro. (Note that there is a strong correlation between the wind time series and the renewable energy time series based on the correlation test.)

| An application of the local joint fluctuation laws-correcting the forecast result of ANN
How to analyze the local feature of the network depends on the specific issue. This paper takes the short-term LMP forecast as an example to demonstrate a typical application of the above conclusions.

Node
In-degree In-degree ratio If there is great confidence in the short-term hydro, wind, and solar power generation forecast, [38][39][40][41][42] it is feasible to forecast the LMP trend according to the above analysis results. For instance, if it is known that the joint fluctuation pattern of the last time is RRFR, given the fluctuation pattern of hydro, wind, and solar power generation of the next step (about 5 minutes later) is RFR, then it can be forecasted that the LMP fluctuation pattern of the next step is most likely to be R, followed by E.
Such conclusions can correct the forecast result of ANN. Figure 16 shows the results of a simple example. A total of 116 consecutive LMPs are forecasted based on ANN, 43 the black solid curves in Figure 16 shows the forecast error. Map {T, R, E, F, B} to {2, 1, 0, −1, −2} to facilitate basic operations. Based on the local law analysis, the evolution trend of LMP is forecasted by three simple methods. Method 1 takes the most likely pattern as the forecast result that is the pattern with the highest in-degree in the local network. And the trend forecast results are shown as the red dotted curves in Figure  16. Method 2 uses the mean calculated from the conversion probability matrix to replace the most likely pattern as the forecast result. And the trend forecast results are shown as the blue dotted curves in Figure 16. Method 3 combines in-CRn of nodes on the basis of Method 2. And the trend forecast results are shown as the green dotted curves in Figure 16. It can be seen that, in general, the correction effect of method 3 is better than that of method 2. In order to observe the consistency of the trend forecast results and the forecast error of ANN, the element-by-element multiplication of their firstorder differences is conducted, and the signs of the results will reflect the consistency. The red "×," blue "+," and green "○"mark the location where the trend forecast results and the forecast error are inconsistent. It is obvious that inconsistent situations are not many. Figure 16 demonstrates that the error of ANN forecast results can be explained and adjusted by the results of this paper.

| CONCLUSIONS AND FUTURE WORK
The joint fluctuation between LMP and renewable energy consumption of ISO-NE are analyzed based on the coarsegrained process and complex networks in this paper.
From a global perspective, the joint fluctuation network follows an exponential distribution and exhibits hierarchical modularity. Through community detection, the joint fluctuation laws are more clearly exposed. And based on the M-HITS algorithm proposed in this paper, the nodes with high statistical analysis value are identified.
From a local perspective, the joint fluctuation network follows a power-law distribution. Different from the preferential attachment rules of scale-free networks, the evolution process of the network has a significant concentration trend. A simple example illustrates the corrective effect of the research results on the forecast results of ANN, which reflects the value of the study.
Although this paper has achieved some results in methods and conclusions, there are still some problems. The research object of this paper is the first-order difference instead of the real value of LMP, hydro, wind, and solar. That is, the level of factors is not considered. In fact, we attempt to screen out the intermediate level data (pattern R/E/F) of LMP for simple analysis and find that the modularity of the joint fluctuation network increased to more than 0.5, and the joint fluctuation laws become more obvious. The change pattern T/B can be used to study unusual fluctuations independently. In addition, the leading and lagging relationship among LMP, hydro, wind, and solar series should be considered. The phase shift of time series is necessary for further analysis to improve the relevant conclusions.