(p124)k分位数

  首先,这个题意思应该是:对一个包含n个元素的集合来说,k分位数是指能把有序集合分成k个等大小集合的k-1个顺序统计量,给出一个能找出某一集合的k分位数的O(nlgk)算法——书上多印了一个“第”字,搞得莫名其妙-_-;

  好吧,既然清楚了,我们来分析以下复杂度,lgk——这说明很有可能要用到分治,再根据递归树模型,我们就可以看出来,大概有这样一个函数:find(int *num,int k,int l,int r,int *result),它的功能是找到第(l+r)/2+1个分隔点对应的顺序统计量,存到result数组里面,然后递归地对左边和右边的区间进行操作。怎么找那个顺序统计量?这一章提到的两种算法都可以,而且这两种算法都顺带将数组分为了小于分割数的一部分和大于分割数的一部分,这保证了这道题分治算法的正确性。

/*
 * source.c
 *
 *  Created on: Feb 16, 2016
 *      Author: wing
 */
#include<stdio.h>
#include<stdlib.h>
int Select(int *num,int m,int l,int r)
{
	int i,j=l-1,tmp;
	if (l==r)
		return num[l];
	for (i=l;i<=r-1;i++)
		if (num[i]<num[r])
		{
			j++;
			tmp=num[j];
			num[j]=num[i];
			num[i]=tmp;
		}
	j++;
	tmp=num[j];
	num[j]=num[r];
	num[r]=tmp;
	if (m-1==j)
		return num[j];
	else
		if (m-1<j)
			return Select(num,m,l,j-1);
		else
			return Select(num,m,j+1,r);
}
int find(int *num,int k,int l,int r,int *result)/*注意哦,这里的k是每个小区间有k个元素^-^*/
{
	result[(l+r)/2]=Select(num,((l+r)/2+1)*k,l*k,(r+2)*k-1);
	if (l==r)
		return 0;
	else
	{
		find(num,k,l,(l+r)/2,result);
		find(num,k,(l+r)/2+1,r,result);
		return 0;
	}
}
int main(void)
{
	int *num,n,k,i,*result;
	scanf("%d",&n);
	scanf("%d",&k);
	num=(int *)malloc(sizeof(int)*n);
	for (i=0;i<n;i++)
		scanf("%d",&num[i]);
	result=(int *)malloc(sizeof(int)*k);
	find(num,n/k,0,k-2,result);
	for (i=0;i<=n/k-2;i++)
		printf("%d ",result[i]);
	return 0;
}


### 分位数格兰杰因果检验概述 分位数格兰杰因果检验是一种扩展的经典格兰杰因果关系检验方法,它允许研究者在不同置下评估因果效应是否存在差异。这种方法特别适用于非对称冲击的时间序列析场景,在金融、经济等领域有广泛应用。 #### 方法论基础 传统的格兰杰因果检验假设变量间的因果关系在整个条件均值范围内是一致的。然而,这种假设可能并不总是成立。分位数格兰杰因果检验则通过考察特定点上的因果关系,揭示更细致的影响模式[^4]。具体而言,该方法利用分位数回归技术来估计滞后项系数,并测试这些系数是否显著异于零。 #### 数学表达形式 对于两个时间序列 \( X_t \) 和 \( Y_t \),如果希望验证 \( X_t \) 是否是 \( Y_t \) 的原因,则可以构建如下模型: \[ Y_t = \alpha_0 + \sum_{i=1}^{p} \beta_i Y_{t-i} + \sum_{j=1}^{q} \gamma_j X_{t-j} + u_t, \] 其中 \( p \) 是自回归阶数,\( q \) 表示外生变量的最大滞后期数。传统格兰杰因果检验关注的是整个条件期望下的参数估计;而分位数版本则是针对给定分位数水平 \( \tau \in (0, 1) \): \[ Q_\tau(Y_t | I_{t-1}) = \alpha_0^\tau + \sum_{i=1}^{p} \beta_i^\tau Q_\tau(Y_{t-i}|I_{t-1}) + \sum_{j=1}^{q} \gamma_j^\tau Q_\tau(X_{t-j}|I_{t-1}), \] 这里 \( Q_\tau(Z|W) \) 表示随机变量 Z 在 W 条件下的第 τ 分位数。通过对各分位数处的 γ 参数进行联合显著性检验即可得出结论。 --- ### Python 实现方案 以下是基于 `statsmodels` 库的一个简单实现框架,展示如何执行分位数 Granger 因果检验: ```python import numpy as np import pandas as pd from statsmodels.regression.quantile_regression import QuantReg def quantile_granger_causality(x_series, y_series, lags_x, lags_y, tau_levels=[0.1, 0.5, 0.9]): """ 执行分位数Granger因果检验 :param x_series: 外部变量序列 :param y_series: 被解释变量序列 :param lags_x: 外部变量最大滞后期数 :param lags_y: 自身变量最大滞后期数 :param tau_levels: 需要检测的分位数列表 :return: 各分位数的结果字典 """ results = {} data = pd.DataFrame({ 'y': y_series, **{f'lag_y{i}': y_series.shift(i) for i in range(1, lags_y+1)}, **{f'lag_x{j}': x_series.shift(j) for j in range(1, lags_x+1)} }).dropna() design_matrix = data.drop(columns=['y']) response_vector = data['y'] for tau in tau_levels: model = QuantReg(response_vector, design_matrix) res = model.fit(q=tau) # 提取γ系数及其显著性信息 gamma_coefs = {col: {'coef': coef, 'pvalue': pval} for col, coef, pval in zip(design_matrix.columns[-lags_x:], res.params[-lags_x:], res.pvalues[-lags_x:])} results[tau] = { 'model_summary': res.summary(), 'gamma_significance': gamma_coefs } return results # 示例调用 if __name__ == "__main__": np.random.seed(42) N = 1000 time_index = pd.date_range(start='2000', periods=N, freq='M') series_X = np.cumsum(np.random.normal(size=N)) series_Y = 0.3 * series_X[:-1] + np.cumsum(np.random.normal(size=N)) result_dict = quantile_granger_causality(pd.Series(series_X), pd.Series(series_Y), lags_x=2, lags_y=2) for tau, details in result_dict.items(): print(f"\nTau Level: {tau}") print(details['model_summary']) # 输出模型摘要 ``` 此脚本定义了一个通用函数用于完成指定分位数层次上的因果推断任务,并返回每一分位数条件下外部因素贡献程度的相关统计量。 --- ### R 实现方案 R 中可以通过组合多个包如 `quantreg`, `vars` 或其他专门处理时间序列工具箱来进行类似操作: ```r library(quantreg) library(zoo) granger_quantile_test <- function(y, x, lag.x, lag.y, taus=c(0.1, 0.5, 0.9)){ df <- cbind(y=y, sapply(seq_len(lag.y), function(k){lag(y,k)}), sapply(seq_len(lag.x), function(k){lag(x,k)}) ) %>% na.omit() names(df)[-(1:(lag.y+1))] <- paste0("X_l", seq_len(lag.x)) models <- list() sig_results <- matrix(nrow=length(taus), ncol=lag.x*2, dimnames=list(NULL,c(paste0("Coef_",seq_len(lag.x)),paste0("PVal_",seq_len(lag.x))))) for(i in seq_along(taus)){ rq_model <- rq(y ~ . ,data=df, tau=taus[i]) coeftest_res <- summary(rq_model)$coefficients models[[as.character(taus[i])]] <<- rq_model sig_results[i, ] <- unlist(c(coeftest_res[grep("^X_l",rownames(coeftest_res)),"Estimate"], coeftest_res[grep("^X_l",rownames(coeftest_res)),"Pr(>|t|)"])) } return(list(models=models,sig=sig_results)) } set.seed(123); tserX <- cumsum(rnorm(1e3)); tserY <- 0.6*tserX[-length(tserX)] + rnorm(length(tserX)-1); res_r <- granger_quantile_test(as.zoo(tserY), as.zoo(tserX), lag.x=3, lag.y=2) print(res_r$sig) ``` 上述代码片段展示了另一种编程环境下开展相同工作的途径——即借助 R 的强大建模能力快速获取结果并加以解读。 --- ### 结语 无论是采用 Python 还是 R 编程环境,都可以灵活运用现有资源开发定制化解决方案以满足实际需求中的复杂情况。值得注意的是,在真实世界应用场景里还需要考虑诸如样本大小限制等因素带来的潜在偏差问题[^4]^。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值