
==== Front
MethodsX
MethodsX
MethodsX
2215-0161
Elsevier

S2215-0161(24)00355-8
10.1016/j.mex.2024.102903
102903
Statistic
Spatial clustering based on geographically weighted multivariate generalized gamma regression
Yasin Hasbi hasbiyasin@live.undip.ac.id
ab
Purhadi purhadi@its.ac.id
a⁎
Choiruddin Achmad choiruddin@its.ac.id
a
a Department of Statistics, Institut Teknologi Sepuluh Nopember, Surabaya 60111 Indonesia
b Department of Statistics, Universitas Diponegoro, Semarang 50275 Indonesia
⁎ Corresponding author. purhadi@its.ac.id
10 8 2024
12 2024
10 8 2024
13 10290322 5 2024
9 8 2024
© 2024 The Author(s)
2024
https://creativecommons.org/licenses/by-nc-nd/4.0/ This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
Geographically Weighted Regression (GWR) is one of the local statistical models that can capture the effects of spatial heterogeneity. This model can be used for both univariate and multivariate responses. However, it should be noted that GWR models require the assumption of error normality. To overcome this problem, we propose a GWR model for generalized gamma distributed responses that can capture the phenomenon of some special continuous distributions. The proposed model is known as Geographically Weighted Multivariate Generalized Gamma Regression (GWMGGR). Parameter estimation is performed using the Maximum Likelihood Estimation (MLE) method optimized with the Bernt-Hall-Hall-Haussman (BHHH) algorithm. To determine the significance of the spatial heterogeneity effect, a hypothesis test was conducted using the Maximum Likelihood Ratio Test (MLRT) approach. We made a spatial cluster based on the estimated model parameters for each response using the k-means clustering method to interpret the obtained results. Some highlights of the proposed method are:• A new model for GWR with multivariate generalized gamma distributed responses to overcome the assumption of normally distributed errors.

• Goodness of fit test to test the spatial effects in GWMGGR model.

• Spatial clustering of districts/cities in Central Java based on three dimensions of educational indicators.

Graphical abstract

Image, graphical abstract

Method name

Geographically Weighted Multivariate Generalized Gamma Regression
Keywords

Spatial heterogeneity
GWMGGR
Maximum likelihood ratio test
K-means cluster
Educational indicators
==== Body
pmcSpecifications tableSubject area:	Mathematics and Statistics	
More specific subject area:	Spatial Statistics, Spatial Heterogeneity, Generalized Linear Models	
Name of your method:	Geographically Weighted Multivariate Generalized Gamma Regression	
Name and reference of original method:	Original method: Multivariate Generalized Gamma RegressionReferences:• H. Yasin, Purhadi, A. Choiruddin, Statistical Inferences for Multivariate Generalized Gamma Regression Model, in: Y. Bee Wah, D. Al-Jumeily OBE, M.W. Berry (Eds.), Data Sci. Emerg. Technol., Springer Nature Singapore, Singapore, 2024: pp. 463–476. https://doi.org/10.1007/978-981-97-0293-0_33.

• H. Yasin, Purhadi, A. Choiruddin, Parameter Estimation and the Goodness-of-fit Test for the Multivariate Generalized Gamma Distribution, in: 2023 Int. Conf. Comput. Control. Informatics Its Appl., 2023: pp. 382–387. https://doi.org/10.1109/IC3INA60834.2023.10285742

	
Resource availability:	The educational indicators and its predictor variables spanning from 2017 to 2021, sourced from BPS of Central Java Province (https://jateng.bps.go.id).	

Background

Point-based spatial regression models that accommodate spatial heterogeneity have been developed, known as Geographically Weighted Regression (GWR) models [[1], [2], [3]]. It has been widely used in several scientific and socio-scientific fields, including environmental health [4], landscape ecology [5], soil quality [6], air quality and remote sensing [7], groundwater management [8], disease patterns [9], [10], and urban heat island effect. GWR models have also been developed for multivariate responses [11,12]. However, it should be noted that these models require the normality assumption of, which is often difficult to satisfied. Several studies have been developed to address this issue, such as GW-Beta Regression [13], GW-Negative Binomial Regression [14,15], GW-Gamma Regression [16], and GW-Weibull Regression [17] models. However, each model is only applicable to a particular distribution.

In this paper, we propose a point-based spatial regression model for response with generalized gamma distribution to accurately describe the positive continuous response variable. The generalized gamma distribution is able to provide several advantages as it has flexible parameters offering a more adaptable structure compared to other continuous distributions, including exponential, gamma, chi-square, log-normal, Erlang, Weibull, half-normal, and Rayleigh [[18], [19], [20]]. Therefore, we extend the GWR model that can be used for several continuous distributions, which is called the Geographically Weighted Multivariate Generalized Gamma Regression (GWMGGR) model. This method is a localized form of the previously developed Multivariate Generalized Gamma Regression (MGGR) model [21]. In this model, parameter estimation is performed - locally so that a unique parameter will be obtained for each location.

Since each location has a unique GWMGGR model regression coefficient, we interpret it by grouping similar locations into a spatial cluster [22,23]. Therefore, the ultimate goal of this study is to create a spatial cluster thematic map based on the GWMGGR model parameters using the k-means clustering method [24]. After that, we use the proposed method for modelling three educational indicators of Central Java: the Mean Years of Schooling (MYS) represents the education development success indicator, and the education participation rate indicator is measured by the Gross Enrolment Rate (GER) and School Enrolment Rate (SER).

Method details

The GWMGGR model is an expansion of the Multivariate Generalised Gamma Regression (MGGR) model. It is built on the premise that the distribution of the response variable adheres to the Multivariate Generalised Gamma Distribution (MGGD). Therefore, in this section we introduce the MGGR model before we explain in depth about parameter estimation, hypothesis testing, and accuracy measurement of the GWMGGR model. Furthermore, we also provide the steps of data analysis using our proposed method.

Model specifications

Multivariate generalized gamma regression (MGGR)

The MGGR model is formulated specifically for MGG-distributed response variables. The probability density function (pdf) of MGG distribution is defined in Eq. (1) [25].(1) f(y|Θ)=(τΓ(λ))Kexp(−(y1−δ1θ1)τ−∑k=2K(yk−yk−1−δkθk)τ)(∏k=1Kθk)(y1−δ1θ1)1−τλ∏k=2K(yk−yk−1−δkθk)1−τλ,

where Θ=[λτθTδT]T, y=[y1y1⋯yK]T, Γ(⋅) is the gamma function, λ,τ>0 are the first and second shape parameters, δT=[δ1δ1⋯δk⋯δK], δk∈R is the threshold parameter for the k-th response, θT=[θ1θ2⋯θk⋯θK], θk>0 is the scale parameter of the k-th response, yk>0, k=1,2,⋯,K, y1+δk<yk; k>1, δ1<y1, and f(y|Θ)=0 for others.

The MGGR model uses the scale parameter as its base, with the “log” function as the link function. It means that changes in the predictor variable are directly associated with changes in the scale parameter, while the shape and threshold parameters remain constant for each observation [21]. In the spatial data framework, the response and predictor variable matrix can be written in Eqs. (2) and (3).(2) Y(n×K)=[y1Ty2T⋯yiT⋯ynT]T,

(3) X(n×(p+1))=[x1Tx2T⋯xiT⋯xnT]T,

where yiT=[y1iy2i⋯yki⋯yKi](1×K), xiT=[1x1ix2i⋯xji⋯xpi](1×(p+1)) i=1,2,⋯,n, yi and xi are vectors of response and predictor variables at the i th location, respectively. Yki is the k-th response variable at the i th location. Therefore, the MGGR model can be written as in Eq. (4).(4) E(Yki)=exp(xiTβk),

where βk=[β0kβ1k⋯βjk⋯βpk]T, j=1,2,⋯,p represents the vector of regression coefficient for the k-th response variable, and p is the number of predictor variables.

The Maximum Likelihood Estimation (MLE) approach is used to estimate the parameters of the MGGR model. This estimation is optimized using the Berndt-Hall-Hall-Hausman (BHHH) method. The MGGR model regression parameters are tested simultaneously using the Maximum Likelihood Ratio Test (MLRT) and partially using the Wald test [21].

Geographically weighted multivariate generalized gamma regression (GWMGGR)

This model is a new extension of the Multivariate Generalized Gamma Regression (MGGR) model, which was developed based on the assumption that the distribution of response variables follows the multivariate generalized gamma distribution [25]. This model is a local form of the MGGR model that applies spatial weights based on latitude and longitude coordinates to the parameter estimation process. This model generates a specific parameter estimator to each point. The general form of the GWMGGR model for each location is displayed in Eq. (5).(5) E(Yki)=exp(xiTβki),

where βki=[β0kiβ1ki⋯βjki⋯βpki]T representing the vector of regression coefficient for the k-th response variable at i th location which is assumed to contain a spatial heterogeneity effect. Therefore, the joint pdf of the GWMGGR model with response variables y1i,y2i,⋯,yKi at the i th location is displayed in Eq. (6).(6) f(yi|ΘGWMRi)=(τiΓ(λi))Kexp(−(y1i−δ1iθ1(xi))τi−∑k=2K(yki−y(k−1)i−δkiθk(xi))τi)(∏k=1Kθk(xi))(y1i−δ1iθ1(xi))1−τiλi∏k=2K(yki−y(k−1)i−δkiθk(xi))1−τiλi,

where ΘGWMRi=[β1iTβ2iT⋯βkiT⋯βKiTλiτiδiT]Tis written as the vector of GWMGGR model parameters at the i th location, λi,τi>0 are the first and second shape parameters at the i th location, δiT=[δ1iδ2i⋯δki⋯δKi], δki∈R is the threshold parameter for the k-th response at the i th location, θk(xi) is defined as the scale parameter of the k-th response at the i th location that is directly affected by the i th predictor, yki>0, y1i+δki<yki;k>1, δ1i<y1i, and f(yi|ΘGWMRi)=0for others.

This model is formulated based on the scale parameter and employs the 'log' function as the link function, so based on the theoretical mean of the MGG distribution, it can be explained that θk(xi) is as follows:

Fork=1,E(Y1i)=θ1(xi)Γ(λi+1τi)Γ(λi)+δ1i=exp(xiTβ1i)or

(7) θ1(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβ1i)−δ1i).

For k=2,E(Y2i)=(θ1(xi)+θ2(xi))Γ(λi+1τi)Γ(λi)+(δ1i+δ2i)=exp(xiTβ2i)

θ2(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβ2i)−(δ1i+δ2i))−θ1(xi),or

θ2(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβ2i)−exp(xiTβ1i)−δ1i).

For k=3,E(Y3i)=(θ1(xi)+θ2(xi)+θ3(xi))Γ(λi+1τi)Γ(λi)+(δ1i+δ2i+δ3i)=exp(xiTβ3i)

θ3(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβ3i)−(δ1i+δ2i+δ3i))−(θ1(xi)+θ2(xi)),or

θ3(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβ3i)−exp(xiTβ2i)−δ3i)

and generally, for k=2,3,⋯,K, θk(xi)can be written in Eq.(8).(8) θk(xi)=Γ(λi)Γ(λi+1τi)(exp(xiTβki)−exp(xiTβ(k−1)i)−δki).

Parameter estimation of GWMGGR model

GWMGGR model parameters are estimated separately at each site by including spatial weighting functions. The generated parameters are varied at different locations. This weight matrix is a representation of the location of one observation relative to another. The weight value increases as the distance between two locations decreases, and vice versa [1]. Suppose the i*-th location is the location to be estimated, and ΘGWMRi*=[β1i*Tβ2i*T⋯βki*T⋯βKi*Tλi*τi*δi*T]T represents the vector of GWMGGR parameters at the i*-th location. First, we construct the log-likelihood function of the GWMGGR model at the i*-th location multiplied by the spatial weights and defined as in Eq. (9).(9) L(ΘGWMRi*)=∑i=1nLi(ΘGWMRi*)=∑i=1nwii*logf(yi|ΘGWMRi*),

where wii* is the weighting value of the i th location with respect to the i*-th location, which depends on the kernel function and the bandwidth used. In this study, we use a gaussian kernel function; consequently, the elements of the weighting matrix can be calculated by Eq. (10).(10) wii*=exp(−12(dii*h)2);i=1,2,⋯,n,

(11) dii*=(ui−ui*)2+(vi−vi*)2,

where dii* is the spatial distance between location i and i* that computed using Euclid distance as in Eq. (11), h is the spatial bandwidth, (ui,vi) and (ui*,vi*) are the geographical coordinate of the i th and i*-th observation. Thus, the spatial weight matrix for the i*-th location is defined in Eq. (12).(12) Wi*=diag(w1i*,w2i*,⋯,wii*,⋯,wni*),

The estimator of ΘGWMRi* can be obtained using the MLE method, namely when the first order condition on Li(ΘGWMRi*) against ΘGWMRi* is equal to zero, and the second order condition on Li(ΘGWMRi*) against ΘGWMRi* is a negative definite matrix. Since the derivation of Li(ΘGWMRi*) with respect to all elements of ΘGWMRi* is not a closed form we use the BHHH iteration method to effectively solves this issue. The primary benefit of this approach is its capacity to generate reliable and effective parameter estimations, as well as its capability to manage intricate nonlinear regression models. The BHHH iteration method synergistically integrates a likelihood criterion for model fitting with an iterative strategy that employs a numerical methodology to optimize the criterion. During each iteration, this method recalculates the parameter estimates by utilizing the inverse of the expected Fisher information matrix, which is computed based on the gradient of the likelihood function. The BHHH method utilizes this strategy to effectively address nonlinearity in the model by iteratively improving the parameter estimates, ultimately leading to convergence towards a solution that maximizes the likelihood [[26], [27], [28]]. In practice, we use this method because it does not require the second derivative to develop the Hessian matrix. The Hessian matrix can be approximated by using the sum of the gradient vector of each observation [29,30]. The gradient vector and the approximation Hessian matrix of the GWMGGR parameters are shown in Eqs. (13) and (14).(13) g(ΘGWMRi*)=∑i=1ngi(ΘGWMRi*),

(14) H(ΘGWMRi*)=−∑i=1ngi(ΘGWMRi*)gi(ΘGWMRi*)T,

where gi(ΘGWMRi*) is the gradient vector of the GWMGGR parameter at the i*-th location for the i th observation and defined in Eq. (15).(15) gi(ΘGWMRi*)=[∂Li(ΘGWMRi*)∂β1i*T⋯∂Li(ΘGWMRi*)∂βKi*T∂Li(ΘGWMRi*)∂λi*∂Li(ΘGWMRi*)∂τi*∂Li(ΘGWMRi*)∂δi*T]T.

The estimation process is repeated until n-th locations to obtain the overall GWMGGR parameters or Θ^GWMR=[Θ^GWMR1Θ^GWMR2⋯Θ^GWMRi⋯Θ^GWMRn]T. Thus, the procedures to obtain the estimator of the GWMGGR model parameters using BHHH iteration can be explained in Fig. 1.Fig. 1 Flowchart of parameter estimation of GWMGGR model.

Fig 1

Hypothesis testing of GWMGGR model

Hypothesis testing of the GWMGGR model parameters simultaneously is carried out with the following hypothesis.

H 0: βjki=0, for each i=1,2,...,n, j=1,2,⋯,p, and k=1,2,⋯,K,

H 1: at least one of βjki≠0.

The statistical test for this hypothesis is generated using the MLRT method. Define that ΩGWMR as the set of GWMGGR model parameters under population and ωGWMR as the set of GWMGGR model parameters under H0, and denoted byΩGWMR={β11,β21,⋯,βki,⋯,βKn,λ1,⋯,λn,τ1,⋯,τn,δ1,⋯,δn},and

ωGWMR={βω011,βω021,⋯,βω0ki,⋯,βω0Kn,λω1,⋯,λωn,τω1,⋯,τωn,δω1,⋯,δωn}.

Furthermore, Ω^GWMR and ω^GWMR are estimated by the MLE method using the parameter estimation procedure as described in Fig. 1. The statistical test of this hypothesis is stated in Eq. (16).(16) GGWMR2=2(L(Ω^GWMR)−L(ω^GWMR)),

where L(Ω^GWMR) and L(ω^GWMR) are the log-likelihood values of the GWMGGR model using the estimated parameters under population and under H0. According to Fotheringham [1], in determining the effective number of parameters in the GWR model, the test statistic GGWMR2 asymptotically for large n follows the chi-squared distribution with the number of effective parameters in the GWMGGR model (Ktr(A)) as the degree of freedom, where A is the projection matrix of the GWMGGR model denoted in Eq. (17). Therefore, the decision to reject H0 can be made if GGWMR2>χ(1−α);Ktr(A)2, where α is the significance level.(17) A(n×n)=[a1Ta2T⋯aiT⋯anT]T,

where aiT=[xiT(XTWiX)−1XTWi](1×n).

Then, to test the model parameters partially, we use the hypothesis H0: βjki=0 versus H1: βjki≠0. The statistical test of this hypothesis is stated in Eq. (18)(18) Zjki=β^jkiVar^(β^jki)→n→∞dN(0,1),

where Var^(β^jki) is the main diagonal element corresponding to β^jki of the matrix −[H(Θ^GWMRi)]−1. So, the null hypothesis is rejected if |Zjki|>Zα/2, where Zα/2 is the quantile of the standard normal distribution.

The relevance of geographic impacts is assessed by comparing the global model (MGGR) with the GWMGGR model using the following hypothesis.

H 0: βjki=βjk;∀i=1,2,...,n;j=0,1,2,⋯,p;k=1,2,⋯,K,

H 1: at least one of βjki≠βjk.

The test statistics of this hypothesis is stated as(19) FGWMR=(GMR2/pK)/(GGWMR2/Ktr(A)),

where GMR2 is the value of the test statistic for testing the hypothesis of the MGGR model simultaneously. H0 is rejected if FGWMR>F(1−α);pK,Ktr(A), where F(1−α);pK,Ktr(A) is the (1−α)-th quantile of the F distribution with pK and Ktr(A) degrees of freedom.

Measures of model fits

Three statistics are used in generalized linear model (GLM) to evaluate model fit: coefficient determination (R2), Deviance, and Akaike Information Criterion (AIC). The coefficient of determination is a statistic designed to assess the correctness of a model and determine how well it fits the data. The coefficient of determination quantifies the extent to which the predictor variables exert an impact on the response variable. A higher coefficient of determination indicates a superior model. The coefficient of determination for the generalized linear model can be computed using the likelihood ratio test (LRT) method outlined in Eq. (20).(20) R2=1−exp(2n(L(ω^)−L(Ω^))),

where n is the number of observations, L(Ω^) denote the log-likelihood of full model (restricted model), and L(ω^) denote the log-likelihood of the null model (unrestricted model) [31,32].

Deviance is a statistical measure used to assess how well a statistical model fits the data. It is commonly employed in statistical hypothesis testing. deviance can be calculated from the log-likelihood model using Eq. (2). The best model is selected based on the smallest Deviance [33].(21) Dev=−2L(Ω^).

AIC is calculated based on the MLE method. The model with the lowest AIC is chosen as the best one. AIC for GWR model can be calculated with Eq. (22) [1].(22) AIC=−2lnL(Ω^)+n+tr(A).

Steps of data analysis 1. Testing the distribution of the response variable with the MGG distribution.

2. Estimate the parameters of the MGGR model.

3. Estimate the parameters of the GWMGGR model using the algorithm described in Fig. 1.

4. Test the GWMGGR model parameters simultaneously using test statistics in Eq. (16).

5. Test the GWMGGR model parameters partially using test statistics in Eq. (18).

6. Test the spatial heterogeneity effect using the test statistics in Eq. (19).

7. Mapping the surface of each regression coefficient of the GWMGGR model.

8. Create spatial cluster based on the regression coefficients of GWMGGR model.

9. Visualize spatial cluster in thematic maps.

10. Make an interpretation and conclusion.

In this study, we analysed the three Educational Indicators in Central Java, Indonesia, using the R programming language. To make the proposed method easier to use with different datasets, we have made the data analysis steps outlined in this study publicly available at the following link (RPubs - Spatial Clustering based on GWMGGR Model).

Method validation

Simulation study

Simulation design

To illustrate the performance of the proposed models, we simulated the GWMGGR and MGGR models on 20 simulated data sets with different degrees of spatial heterogeneity on the local parameter surface. We used a spatial layout consisting of a 10×10 grid (vertical vs horizontal coordinates). We set up three simulation scenarios for one predictor based on three levels of spatial heterogeneity, namely zero, low, and high, under the following rules [34]:

Simulation Scenario 1 (zero spatial heterogeneity):

In this scenario all regression coefficients are assumed to be equal for all points, which are as followsβ01i=2.45,β02i=3.51,andβ03i=3.95

β11i=0.0045,β12i=0.0012,andβ13i=0.0014

Simulation Scenario 2 (low spatial heterogeneity):

In scenario 2, the regression coefficients for each location are generated with the following rulesβ01i=2.45,β02i=3.51andβ03i=3.95

β11i=log(1+(112)(ui+vi))/100,β12i=log(1+(112)(ui+vi))/75,andβ13i=log(1+(112)(ui+vi))/85.

Simulation Scenario 3 (high spatial heterogeneity):β01i=2.45,β02i=3.51,andβ03i=3.95

β11i=log(1+(1324)(36−(6−ui2)2)(36−(6−vi2)2))/100,

β12i=log(1+(1324)(36−(6−ui2)2)(36−(6−vi2)2))/75,and

β13i=log(1+(1324)(36−(6−ui2)2)(36−(6−vi2)2))/85,

where vi and ui are the vertical and horizontal coordinate values, which increases by 1 with each increment. The algorithm used to simulate the model is given in Table 1.Table 1 Algorithm for GWMGGR simulation.

Table 1Define : The spatial layout consists of a regular 10×10 grid as geographical coordinates ui=(ui,vi),
so the number of locations is n = 100
Number of response K = 3, and number of predictor p = 1.
Assume : The true GWMGGR regression parameters β1i,β2i,andβ3i(for three response variables) following the three designs to be simulated.
The MGG distribution parameters λ=1.25,τ=4.5,andδT=[3,10,8].
Distribution of predictor variables X1∼U(5,30)
Data Generating Process
fori=1,2,⋯,ndo
generatexi=[1x1i]T
computeθ1(xi) using Eq. (7)
θ2(xi) and θ3(xi) using Eq. (8)
generateyi∼MGG(λ,τ,θi,δ), where θi=[θ1(xi)θ2(xi)θ3(xi)]T
end for
Apply GWMGGR to yon x1	

Simulation result

We present some statistics of the MGGR and GWMGGR models with three spatial heterogeneity scenarios (zero, low, and high) in Table 2. There are several points that can be concluded:1. The optimum bandwidth of the GWMGGR model will be smaller for data with higher levels of spatial heterogeneity. The smaller the bandwidth, each location will only be influenced by several neighbouring locations. Thus, the GWMGGR model parameters will have a higher variation.

2. The smaller bandwidth used, the effective parameter will increase and the AIC of the GWMGGR model also tends to become larger.

3. The high level of spatial heterogeneity causes the log-likelihood of the MGGR model to be smaller, which leads to the MGGR model being inaccurate. Conversely, the accuracy of the GWMGGR model will be improved. This can be seen in the coefficient of determination of spatial heterogeneity, the decrease in Deviance, and the F-Stat test statistics that are progressively increasing.

4. Probability (p-value) of the model similarity test (MGGR and GWMGGR) will become smaller as the level of spatial heterogeneity increases. This means that it can be inferred that the performance of the GWMGGR model is superior to be able to overcome the presence of spatial heterogeneity.

Table 2 Simulation results.

Table 2Average of the Statistics	Level of Spatial Heterogeneity	
Zero	Low	High	
Optimum Bandwidth	12.661	2.400	1.623	
n-eff	6.606	27.963	48.944	
AIC MGGR	1370.100	1528.168	1678.617	
AIC GWMGGR	1463.287	1574.979	1644.743	
Log-likelihood MGGR	−679.050	−758.084	−833.308	
Log-likelihood GWMGGR	−678.341	−723.508	−747.899	
Determination of Spatial Heterogeneity (%)	0.628	26.110	52.798	
Deviance MGGR	1358.100	1516.168	1666.617	
Deviance GWMGGR	1356.682	1447.016	1495.799	
Decrease in Deviance	1.418	69.151	170.818	
F-Stat GWMGGR	2.035	8.474	13.045	
p-value F-Test	2.030×10−1	1.963×10−3	1.993×10−4	

Application on three educational indicators data

Data set

We use secondary data sourced from BPS Central Java Province, which consists of 35 districts/cities for the period 2017 to 2021. The response variable consists of three education development indicators Mean Years of Schooling (MYS), Gross Enrolment Rate (GER) and School Enrolment Rate (SER) which are multivariate generalized gamma distributed [25]. The predictor variables that are expected to influenced the three indicators are Income per capita, Percentage of Households Accessing Proper Sanitation, and Student-Teacher Ratio [35]. Table 3 shows a summary of the research data. As an overview, the average MYS in Central Java Province in that period reached 7.77 years, which means that the education of the population aged 25 years and over only reached the level of first grade of junior high school or dropped out of school in second grade. This is even though only 17% of districts and cities have achieved the 9-year compulsory education target. The average SER for 16–18 years old is 71.50, which means that only around 71.50% of the population aged 16–18 years old are currently in school, and this still does not match the national target. The average GER at the senior high school level is 85.73, which means that only around 85.73% of the population aged 16–18 are currently attending senior high school, which again falls short of the national target of 88.39%.Table 3 Descriptive statistics of research data.

Table 3Variable	Description	Mean	SD	Min	Max	
Response	MYS (Y1) (years)	7.773	1.212	6.180	10.900	
	SER for 16–18 years (Y2) (%)	71.501	9.011	49.560	91.390	
	GER at the high school (Y3) (%)	85.734	14.941	52.980	121.910	
Predictor	GRDP per capita (X1) (million Rupiah)	28.386	17.966	12.370	87.360	
	Percentage of Households Accessing Proper Sanitation (X2) (%)	77.954	17.292	9.240	98.070	
	Student-Teacher Ratio (X3)	18.063	2.934	13.000	39.000	

Modelling educational indicators using MGGR

We first tested the distribution of response variables simultaneously using the Kolmogorov-Smirnov test. Based on the three response variables in this study, the test statistic was 0.0597 with a p-value of 0.5616 (more than 0.05). Thus, it can be concluded that the response variable data follows the MGG distribution [25]. Furthermore, we estimate the parameters of the MGGR model using the likelihood method and the BHHH algorithm. The results of the MGGR model parameter estimation using three predictor variables are presented in Table 5. Using a significant level of 5%, it can be concluded that all predictors significantly affect all responses. The statistical significance of the model can be evaluated simultaneously by utilizing Wilk's likelihood ratio statistics obtained from the MLRT (refer to Table 4). The test statistic is estimated at 122.0063, and the p-value is 0.0000. It implies that all three factors have a significant influence on the response variables collectively. Thus, mathematically, the MGGR model can be presented in Eq. (23).(23) E(y1)^=exp(1.8108+0.0044x1+0.0028x2−0.0059x3)E(y2)^=exp(4.1557+0.0011x1+0.0027x2−0.0095x3)E(y3)^=exp(4.1069+0.0013x1+0.0047x2−0.0058x3)

Table 4 Test of MGGR parameters simultaneously.

Table 4Source	Deviance	GMR2	df	χ(0.95);92	p-value	
Null Model	3158.794	122.0063	9	16.919	0.0000	
Residuals	3036.787					

Table 5 Regression coefficients of MGGR model.

Table 5Par.	MYS (Y1)	SER (Y2)	GER (Y3)	
Est.	Z	p-value	Est.	Z	p-value	Est.	Z	p-value	
β0	1.8108	1.19E+10	0.000*	4.1557	9.33E+18	0.000*	4.1069	3.69E+16	0.000*	
β1	0.0044	1.81E+06	0.000*	0.0011	1.54E+14	0.000*	0.0013	7.34E+11	0.000*	
β2	0.0028	1.58E+06	0.000*	0.0027	4.73E+14	0.000*	0.0047	3.61E+12	0.000*	
β3	−0.0059	−2.41E+06	0.000*	−0.0095	−1.32E+15	0.000*	−0.0058	−3.27E+12	0.000*	
⁎ Significant at α=5%.

Modelling educational indicators using GWMGGR

This model is heavily impacted by the spatial weight matrix, whose value depends on the distance between the observed coordinates. The choice of the kernel function and bandwidth significantly affects the most suitable weight matrix. Hence, we employed the grid-search technique to determine the optimal bandwidth within the range of minimum and maximum distances between places, as illustrated in Table 6. The optimal bandwidth was selected using the criteria of minimising Deviance and AIC, and maximising R2. Consequently, the optimal bandwidth for this case study was determined to be 0.25. By using the algorithm in Fig. 1, a summary of GWMGGR parameter estimates is obtained in Table 7 (the complete coefficients are given in Appendix 1). To see the influence spread of each predictor on each response, it is shown in the surface map in Fig. 2. The inference of the GWMGGR model was conducted with three hypothesis tests. The first hypothesis is about the simultaneous test of model parameters. Based on Table 8, the test statistic GGWMR2 = 124.049 and p-value = 0.0051, it can be inferred that the three predictor variables in the GWMGGR model simultaneously have a significant effect on the three educational indicators. The second hypothesis is about partial model parameter test. Because the GWMGGR model produces unique parameter estimates for each location, hypothesis testing is carried out for each parameter at each location using the test statistics in Eq. (18). For example, the partial hypothesis test for Wonosobo District (the 7th location) is presented in Table 9. According to this table, it can be inferred that all predictors significantly affect the three educational indicators in Wonosobo.Table 6 Bandwidth selection in GWMGGR using grid-search.

Table 6bandwidth	R2 (%)	Loglikelihood	Deviance	AIC	
0.20	89.229	−1384.417	2768.834	3054.845	
0.25	90.148*	−1376.621*	2753.242*	3014.788*	
0.30	84.297	−1417.406	2834.813	3079.944	
0.35	85.765	−1408.821	2817.643	3051.139	
0.40	81.187	−1433.217	2866.435	3091.344	
0.45	74.738	−1459.009	2918.018	3136.405	
0.50	71.117	−1470.729	2941.459	3154.796	
⁎ Optimum bandwidth.

Table 7 Summary of regression coefficients in GWMGGR model.

Table 7Response	Parameter	Minimum	Q1	Median	Mean	Q3	Maximum	
MYS (Y1)	β01	1.477809	1.703308	1.742911	1.742306	1.810728	1.940631	
β11	−0.004687	0.004410	0.005448	0.005529	0.006561	0.014184	
β21	−0.000486	0.001837	0.002266	0.002141	0.002764	0.003334	
β31	−0.010486	−0.006172	−0.000819	−0.001296	0.001096	0.013444	
SER (Y2)	β02	3.631474	4.116476	4.155845	4.118853	4.181724	4.350242	
β12	−0.000622	0.000932	0.001476	0.002423	0.002907	0.009549	
β22	−0.001331	0.000978	0.002351	0.002125	0.002969	0.007367	
β32	−0.021871	−0.009521	−0.008310	−0.006533	−0.005683	0.019608	
GER (Y3)	β03	3.624337	4.106884	4.134112	4.202760	4.278730	4.917275	
β13	−0.001555	0.000677	0.001282	0.002528	0.003388	0.012819	
β23	−0.003494	0.002249	0.004877	0.003697	0.005473	0.007426	
β33	−0.022958	−0.013823	−0.008134	−0.008315	−0.005828	0.019882	

Table 8 Test of GWMGGR parameters simultaneously.

Table 8Source	Deviance	GGWMR2	df	χ(0.95);86.5462	p-value	
Null Model	3036.787	124.049	86.546	109.263	0.0051	
Residuals	2753.242					

Table 9 Regression coefficients in Wonosobo District (the 7th location).

Table 9Par	MYS (Y1)	SER (Y2)	GER (Y3)	
Est.	Z	p-value	Est.	Z	p-value	Est.	Z	p-value	
β0	1.8108	1.19E+10	0.000*	4.1557	8.36E+17	0.000*	4.1069	3.59E+16	0.000*	
β1	0.0044	1.81E+06	0.000*	0.0012	1.54E+13	0.000*	0.0013	7.06E+11	0.000*	
β2	0.0028	1.61E+06	0.000*	0.0027	4.96E+13	0.000*	0.0046	3.46E+12	0.000*	
β3	−0.0059	−2.42E+06	0.000*	−0.0094	−1.18E+14	0.000*	−0.0059	−3.23E+12	0.000*	
⁎ Significant at α=5%.

Fig. 2 Surface of the GWMGGR Coefficients.

Fig 2

The third hypothesis is to test the suitability of the MGGR and GWMGGR models to determine whether there is an effect of spatial heterogeneity. The test statistics are calculated using Eq. (19), and the results are presented in Table 10. The table inform us that FGWMR = 9.4579 and p-value = 0.0000, which indicates that H0 is rejected. Thus, it can be said that the GWMGGR model is more suitable than the MGGR model, or there is a significant effect of spatial heterogeneity on modelling educational indicators in Central Java. This decision is also supported by the higher coefficient of determination and the lower AIC of the GWMGGR model compared to the MGGR model, as shown in Table 11. The GWMGGR model can explain the variability up to 90.148% compared to the MGGR model which can only explain it by 50.201%. It is concluded that the GWMGGR model is effective in overcoming the problem of spatial heterogeneity.Table 10 Test of model similarity (GWMGGR vs MGGR).

Table 10Model	G2	df	F	F-table	p-value	
MGGR	122.0063	9	9.4579	1.9899	0.0000	
GWMGGR	124.0490	86.5462				

Table 11 Model Comparison.

Table 11Model	R2 (%)	AIC	
GWMGGR	90.148*	3014.788*	
MGGR	50.201	3060.787	
⁎ Best model.

Spatial clustering based on the regression coefficient of GWMGGR model

After obtaining the best parameters of GWMGGR, further model interpretation can be done by creating a spatial cluster of districts/cities based on the regression coefficient of the GWMGGR model at each location. We used the k-means clustering method to form spatial clusters, and to determine the optimal number of clusters, we applied the silhouette method [36]. Based on Fig. 3, the best number of clusters is 3 with a silhouette value of 0.3403. The cluster membership of each cluster in detail are presented in Table 12. Visually, the spatial cluster of the GWMGGR model parameters is presented in a thematic map in Fig. 4. The map shows that districts/cities in Central Java are grouped into three clusters, each with unique characteristics. These results indicate that to improve the achievements of these three educational indicators, a special effort is needed in accordance with the characteristics of each district/city clusters.Fig. 3 Silhouette score based on the number of clusters.

Fig 3

Table 12 Spatial cluster membership based on regression coefficients of GWMGGR model.

Table 12Cluster	Members	District/City	
1	4	Tegal, Pemalang, Pekalongan City, and Tegal City	
2	25	Banyumas, Purbalingga, Banjarnegara, Kebumen, Purworejo, Wonosobo, Magelang, Boyolali, Klaten, Sukoharjo, Wonogiri, Karanganyar, Sragen, Grobogan, Jepara, Demak, Semarang, Temanggung, Kendal, Batang, Pekalongan, Magelang City, Surakarta City, Salatiga City, and Semarang City	
3	6	Cilacap, Brebes, Blora, Rembang, Pati, and Kudus	

Fig. 4 Spatial cluster of district/city in Central Java.

Fig 4

Conclusions

This study is mostly about creating point-based spatial models for regression modelling using multivariate generalized gamma responses to deal with data that isn't normal. We introduce the GWMGGR model, which is a local variation of the MGGR model. The proposed model uses the geographical coordinates of observation locations as a weighting factor to represent spatial heterogeneity. This model will generate different parameter estimates for each location, allowing us to interpret the model as a spatial cluster map to group locations with similar characteristics. The results on modelling educational indicator data in Central Java concludes that the GWMGGR model can overcome spatial heterogeneity, as evidenced by its increased accuracy compared to the MGGR model. The parameters of the GWMGGR model categorize districts and cities in Central Java into three groups, each necessitating specific measures to enhance the achievement of educational indicators. The findings of this study can inform further research on the MGGR model, with a particular focus on the presence of temporal and spatio-temporal heterogeneity factors.

Limitations

None.

Ethics statements

The data used in this research are secondary data derived from the official website of BPS of Central Java Provinces Indonesia (https://jateng.bps.go.id/).

CRediT authorship contribution statement

Hasbi Yasin: Conceptualization, Methodology, Software, Writing – original draft. Purhadi: Conceptualization, Methodology, Writing – review & editing, Supervision. Achmad Choiruddin: Methodology, Writing – review & editing, Validation, Supervision.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A

Tables A1, and A2Tabel A1. Regression coefficient of GWMGGR model.

Tabel A1No	MYS (Y1)	SER (Y2)	GER (Y3)	
Intercept	X1	X2	X3	Intercept	X1	X2	X3	Intercept	X1	X2	X3	
1	1.8406	−0.0018	0.0028	−0.0017	4.0863	0.0019	0.0020	−0.0077	4.4250	0.0011	−0.0025	0.0084	
2	1.7273	0.0044	0.0022	−0.0002	4.1853	0.0011	0.0027	−0.0117	4.1365	0.0016	0.0049	−0.0060	
3	1.8563	0.0021	0.0005	0.0015	4.0059	0.0010	0.0074	−0.0219	4.0327	0.0045	0.0074	−0.0137	
4	1.8108	0.0044	0.0028	−0.0059	4.1557	0.0011	0.0029	−0.0094	4.1069	0.0013	0.0045	−0.0059	
5	1.8107	0.0048	0.0030	−0.0064	4.1557	0.0015	0.0035	−0.0096	4.1069	0.0012	0.0052	−0.0058	
6	1.8108	0.0049	0.0031	−0.0066	4.1557	0.0019	0.0034	−0.0094	4.1069	0.0007	0.0052	−0.0062	
7	1.8108	0.0044	0.0028	−0.0059	4.1557	0.0012	0.0029	−0.0094	4.1069	0.0013	0.0046	−0.0059	
8	1.8079	0.0054	0.0033	−0.0072	4.1511	0.0024	0.0034	−0.0119	4.1141	0.0021	0.0058	−0.0134	
9	1.4778	0.0070	0.0023	0.0130	4.1707	0.0012	0.0024	−0.0065	4.3339	0.0014	0.0049	−0.0185	
10	1.5056	0.0066	0.0024	0.0116	4.1565	0.0008	0.0026	−0.0057	4.2631	0.0005	0.0054	−0.0149	
11	1.5530	0.0069	0.0027	0.0082	4.1928	0.0025	0.0023	−0.0083	4.1126	0.0027	0.0050	−0.0071	
12	1.7276	0.0036	0.0027	0.0002	4.0903	0.0027	0.0030	−0.0081	4.1112	0.0008	0.0066	−0.0152	
13	1.4891	0.0068	0.0018	0.0134	4.1427	−0.0006	0.0029	−0.0043	4.2943	−0.0014	0.0055	−0.0151	
14	1.7429	0.0059	0.0022	−0.0021	4.1740	−0.0006	0.0030	−0.0063	4.1941	−0.0013	0.0057	−0.0092	
15	1.7201	0.0054	0.0021	−0.0007	4.2126	0.0012	0.0032	−0.0133	4.1589	0.0007	0.0067	−0.0140	
16	1.7589	0.0028	0.0008	0.0032	4.3502	0.0003	0.0003	−0.0090	4.6737	0.0010	0.0013	−0.0229	
17	1.9387	0.0027	−0.0005	0.0004	4.1450	0.0005	0.0010	−0.0017	4.5172	0.0005	0.0011	−0.0111	
18	1.9406	−0.0047	0.0017	−0.0007	4.1846	0.0050	0.0022	−0.0125	4.5320	−0.0016	0.0023	−0.0143	
19	1.8558	0.0029	0.0019	−0.0090	4.2058	0.0064	−0.0013	−0.0012	4.9173	0.0069	−0.0035	−0.0230	
20	1.7835	0.0048	0.0006	0.0056	4.0161	0.0009	0.0022	0.0023	4.1440	0.0003	0.0046	0.0007	
21	1.7254	0.0055	0.0030	−0.0015	4.0011	0.0037	0.0010	0.0030	4.0461	0.0020	0.0043	0.0021	
22	1.7084	0.0062	0.0022	−0.0008	4.1722	0.0007	0.0035	−0.0123	4.1301	0.0003	0.0068	−0.0132	
23	1.6982	0.0064	0.0023	−0.0004	4.1789	0.0014	0.0023	−0.0083	4.1341	0.0013	0.0052	−0.0088	
24	1.7970	0.0049	0.0027	−0.0055	4.1589	0.0030	0.0004	−0.0057	4.1026	0.0041	0.0042	−0.0078	
25	1.7809	0.0066	0.0023	−0.0059	4.1662	0.0048	0.0005	−0.0093	4.0927	0.0117	0.0031	−0.0148	
26	1.8088	0.0058	0.0028	−0.0072	4.1557	0.0024	0.0009	−0.0070	4.1068	0.0067	0.0024	−0.0077	
27	1.7277	0.0142	0.0020	−0.0098	3.6315	0.0076	−0.0002	0.0196	3.6243	0.0099	0.0011	0.0199	
28	1.6526	0.0124	0.0014	−0.0030	4.2261	0.0095	−0.0009	−0.0109	4.2517	0.0051	0.0009	−0.0065	
29	1.8581	0.0065	0.0014	−0.0091	3.9868	0.0020	0.0033	−0.0091	4.4726	0.0023	−0.0018	−0.0049	
30	1.8107	0.0050	0.0031	−0.0067	4.1558	0.0028	0.0025	−0.0093	4.1066	0.0026	0.0049	−0.0097	
31	1.6747	0.0061	0.0024	0.0016	4.1898	−0.0001	0.0027	−0.0070	4.1665	−0.0003	0.0055	−0.0081	
32	1.7119	0.0050	0.0025	0.0000	4.1753	0.0009	0.0038	−0.0128	4.1417	0.0007	0.0066	−0.0136	
33	1.6868	0.0062	0.0022	0.0007	4.2241	0.0013	0.0015	−0.0079	4.0932	0.0007	0.0056	−0.0045	
34	1.7368	0.0142	0.0023	−0.0105	3.8660	0.0062	0.0001	0.0070	3.9203	0.0128	0.0022	−0.0032	
35	1.6337	0.0091	0.0010	0.0019	3.7789	0.0064	0.0009	0.0070	4.3189	0.0042	−0.0022	0.0026	

Table A2. Description of location number.

Table A2No	District/City	No	District/City	
1	Cilacap	19	Kudus	
2	Banyumas	20	Jepara	
3	Purbalingga	21	Demak	
4	Banjarnegara	22	Semarang	
5	Kebumen	23	Temanggung	
6	Purworejo	24	Kendal	
7	Wonosobo	25	Batang	
8	Magelang	26	Pekalongan	
9	Boyolali	27	Pemalang	
10	Klaten	28	Tegal	
11	Sukoharjo	29	Brebes	
12	Wonogiri	30	Magelang	
13	Karanganyar	31	Surakarta City	
14	Sragen	32	Salatiga City	
15	Grobogan	33	Semarang City	
16	Blora	34	Pekalongan City	
17	Rembang	35	Tegal City	
18	Pati			

Appendix B Supplementary materials

Image, application 1

Data availability

Data will be made available on request.

Acknowledgments

The first author would like to express deep appreciation for the financial support provided by Balai Pembiayaan Pendidikan Tinggi (BPPT) or Central for Higher Education Funding and Lembaga Pengelola Dana Pendidikan (LPDP) under the Ministry of Education, Culture, Research, and Technology of the Republic of Indonesia.

Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.mex.2024.102903.
==== Refs
References

1 Fotheringham A.S. Brunsdon C. Charlton M. Geographically Weighted Regression: The Analysis of Spatially Varying Relationships 2002 John Wiley & Sons Inc. Chichester
2 De La Hoz-M J. Mendes S. Fernandez-Gómez M.J. GeoWeightedModel : an R-shiny package for geographically weighted models SoftwareX 20 2022 101250 10.1016/j.softx.2022.101250
3 Lu B. Hu Y. Yang D. Liu Y. Liao L. Yin Z. Xia T. Dong Z. Harris P. Brunsdon C. Comber L. Dong G. GWmodelS: a software for geographically weighted models SoftwareX 21 2023 101291 10.1016/j.softx.2022.101291
4 Yoneoka D. Saito E. Nakaoka S. New algorithm for constructing area-based index with geographical heterogeneities and variable selection : an application to gastric cancer screening Sci. Rep. 6 2016 1 7 10.1038/srep26582 28442746
5 Yu T. Bao A. Xu W. Guo H. Jiang L. Zheng G. Yuan Y. NZABARINDA V. Exploring variability in landscape ecological risk and quantifying its driving factors in the amu darya delta Int. J. Environ. Res. Public Health 17 2020 10.3390/ijerph17010079
6 Zhang C. Tang Y. Xu X. Kiely G. Towards spatial geochemical modelling: use of geographically weighted regression for mapping soil organic carbon contents in Ireland Appl. Geochem. 26 2011 1239 1248 10.1016/j.apgeochem.2011.04.014
7 Wang Q. Feng H. Feng H. Yu Y. Li J. Ning E. The impacts of road traffic on urban air quality in Jinan based GWR and remote sensing Sci. Rep. 11 2021 1 9 10.1038/s41598-021-94159-8 33414495
8 Koh E.-H. Lee E. Lee K.-K. Application of geographically weighted regression models to predict spatial characteristics of nitrate contamination: implications for an effective groundwater management strategy J. Environ. Manage. 268 2020 110646 10.1016/j.jenvman.2020.110646
9 Isazade V. Qasimi A.B. Dong P. Kaplan G. Isazade E. Integration of Moran's I, geographically weighted regression (GWR), and ordinary least square (OLS) models in spatiotemporal modeling of COVID-19 outbreak in Qom and Mazandaran Provinces Iran, Model. Earth Syst. Environ. 9 2023 3923 3937 10.1007/s40808-023-01729-y
10 Brunton L.A. Alexander N. Wint W. Ashton A. Broughan J.M. Using geographically weighted regression to explore the spatially heterogeneous spread of bovine tuberculosis in England and Wales Stoch. Environ. Res. Risk Assess. 31 2017 339 352 10.1007/s00477-016-1320-9
11 Chen V.Y.-J. Yang T.-C. Jian H.-L. Geographically weighted regression modeling for multiple outcomes Ann. Am. Assoc. Geogr. 112 2022 1278 1295 10.1080/24694452.2021.1985955
12 Harini S. Purhadi Parameter estimation of Multivariate Geographically Weighted Regression model using matrix laboratory 2012 Int. Conf. Stat. Sci. Bus. Eng. 2012 1 4 10.1109/ICSSBE.2012.6396622
13 da Silva A.R. de Oliveira Lima A. Geographically weighted beta regression Spat. Stat. 21 2017 279 303 10.1016/j.spasta.2017.07.011
14 Ricardo Da Silva A. Carvalho T. Rodrigues V. Da Silva A.R. Rodrigues T.C.V Geographically weighted negative binomial regression-incorporating overdispersion Stat. Comput. 24 2014 769 783 10.1007/s11222-013-9401-9
15 Yasin H. Suryani I. Kartikasari P. Graphical interface of geographically weighted negative binomial regression (GWNBR) model using R-Shiny J. Phys. Conf. Ser. 2021 1943 10.1088/1742-6596/1943/1/012155
16 Yasin H. Inayati S. Setiawan 3-Parameter gamma regression model for analyzing human development index of central Java Province BAREKENG J. Ilmu Mat. Dan Terap. 16 2022 171 180 10.30598/barekengvol16iss1pp171-180
17 Suyitno N.W.W.Sari Parameter estimation of mixed geographically weighted Weibull regression model J. Phys. Conf. Ser. 1277 2019 1 11 10.1088/1742-6596/1277/1/012046
18 Sanchez R. Mackenzie S.A. Information thermodynamics of cytosine DNA methylation PLoS ONE 11 2016 1 20 10.1371/journal.pone.0150427
19 Shanker R. Shukla K.K. On modeling of lifetime data using three-parameter generalized lindley and generalized gamma distributions Biom. Biostat. Int. J. 4 2016 283 288 10.15406/bbij.2016.04.00117
20 Diantini N.L.S. Purhadi A.C. Parameter estimation and hypothesis testing on three parameters log normal regression AIP Conf. Proc 2023 AIP Publishing LLC 30024 10.1063/5.0104443
21 Yasin H. Purhadi A.C. Statistical inferences for multivariate generalized gamma regression model Bee Wah Y. Al-Jumeily OBE D. Berry M.W. Data Sci. Emerg. Technol 2024 Springer Nature Singapore Singapore 463 476 10.1007/978-981-97-0293-0_33
22 Matthews S.A. Yang T.-C. Mapping the results of local statistics: using geographically weighted regression Demogr. Res. 26 2012 151 166 10.4054/DemRes.2012.26.6 25578024
23 Mennis J. Mapping the results of geographically weighted regression Cartogr. J. 43 2006 171 179 10.1179/000870406x114658
24 Sumarah J. Wulandari A.T. Analysis of the K-means algorithm for clustering school participation rates in Central Java KnE Soc. Sci. 8 2023 10.18502/kss.v8i9.13320
25 Yasin H. Purhadi A.C. Parameter estimation and the goodness-of-fit test for the multivariate generalized gamma distribution 2023 Int. Conf. Comput. Control. Informatics Its Appl. 2023 382 387 10.1109/IC3INA60834.2023.10285742
26 Wooldridge J.M. Econometric Analysis of Cross Section and Panel Data 2010 The MIT Press London
27 Greene W.H. Econometric Analysis 7th ed. 2012 Pearson Education New Jersey
28 Berndt E.K. Hall B.H. Hall R.E. Hausman J.A. Estimation and inference in nonlinear structural models Ann. Econ. Soc. Meas. 3 1974 653 665
29 Purhadi A.R. Wenur G.H. Geographically weighted three-parameters bivariate gamma regression and its application Symmetry 13 2021 1 17 10.3390/sym13020197
30 Wenur G.H. Purhadi A.S. Three-parameter bivariate gamma regression model for analyzing under-five mortality rate and maternal mortality rate J. Phys. Conf. Ser. 1538 2020 1 10 10.1088/1742-6596/1538/1/012054
31 Magee L. R2 measures based on wald and likelihood ratio joint significance tests Am. Stat. 44 1990 250 253 10.2307/2685352
32 Zhang D. A coefficient of determination for generalized linear models Am. Stat. 71 2017 310 316 10.1080/00031305.2016.1256839
33 Collett D. Modelling Survival Data in Medical Research 4th ed. 2023 Chapman and Hall/CRC New York 10.1201/9781003282525
34 Fotheringham A.S. Yang W. Kang W. Multiscale geographically weighted regression (MGWR) Ann. Am. Assoc. Geogr. 107 2017 1247 1265 10.1080/24694452.2017.1352480
35 Gumus S. Kayhan S. The relationship between economic growth and school enrollment rates: time series evidence from Turkey Educ. Policy Anal. Strateg. Res. 7 2012 24 38 https://eric.ed.gov/?id=EJ1127574
36 Kaufman L. Rousseeuw P.J. Finding Groups in Data: An Introduction to Cluster Analysis 2005 John Wiley & Sons Inc. New York 10.1002/9780470316801
