diff options
Diffstat (limited to 'R/interpweights.R')
| -rw-r--r-- | R/interpweights.R | 36 |
1 files changed, 16 insertions, 20 deletions
diff --git a/R/interpweights.R b/R/interpweights.R index 1e9d50c5..f4f6cb50 100644 --- a/R/interpweights.R +++ b/R/interpweights.R @@ -1,15 +1,17 @@ interpweights <- function(w, v1, v2){
#Given L=(w,v1), compute neww such that newL=(new,v2)=L in distribution
cumw <- cumsum(w)
- neww <- splinefun(v1,cumw,method= "monoH.FC")(v2,deriv=1)
+ neww <- splinefun(v1, cumw, method= "monoH.FC")(v2, deriv=1)
#neww <- diff(newcumw)
interpweights <- neww/sum(neww)
return(interpweights)
}
interpvalues <- function(w, v, neww){
+ ## Given a distribution D=(w,v), compute new weights
+ ## such that Dnew=(neww, newv) equals D in distribution
cumw <- cumsum(w)
- test <- splinefun(v, cumw, method="monoH.FC")
+ cdf <- splinefun(v, cumw, method="monoH.FC")
eps <- 1e-3
newv <- rep(0, length(neww))
cumneww <- cumsum(neww)
@@ -17,20 +19,20 @@ interpvalues <- function(w, v, neww){ for(i in 1:length(neww)){
iter <- 0
## do binary search
- hi <- test(1)
+ hi <- cdf(1)
lo <- mid
if(hi < cumneww[i]){
newv[i] <- hi
next
}
- if(test(lo) > cumneww[i]){
+ if(cdf(lo) > cumneww[i]){
newv[i] <- lo
next
}
mid <- (lo+hi)/2
iter <- 0
- while(abs(test(mid)-cumneww[i])>eps){
- if(test(mid)>cumneww[i]){
+ while(abs(cdf(mid)-cumneww[i])>eps){
+ if(cdf(mid)>cumneww[i]){
hi <- mid
}else{
lo <- mid
@@ -42,6 +44,14 @@ interpvalues <- function(w, v, neww){ return(newv)
}
+interpvalues.distr <- function(w, v, neww){
+ ## same as interpvalues, but using the distr
+ ## package. need to check how good it is
+ require(distr)
+ D <- DiscreteDistribution(w, v)
+ return( q(D)(neww) )
+}
+
adjust_scenario <- function(scenario, epsilon){
1-(1-scenario)^(1/(1+epsilon))
}
@@ -113,17 +123,3 @@ inverse <- function(f,Df,x, x0=x){ }
return(x1)
}
-setwd("W:/corpCDOS/NextGen")
-lcdx10 <- read.table("lcdx10",sep=",",header=T)
-lcdx12 <- read.table("lcdx12",sep=",",header=T)
-hy10_5y <- read.table("hy10_5y",sep=",",header=T)
-hy10_7y <- read.table("hy10_7y",sep=",",header=T)
-
-lcdx10_pv <- read.table("lcdx10_pv",sep=",",header=T)
-lcdx12_pv <- read.table("lcdx12_pv",sep=",",header=T)
-hy10_5y_pv <- read.table("hy10_5y_pv",sep=",",header=T)
-hy10_7y_pv <- read.table("hy10_7y_pv",sep=",",header=T)
-clipw(interpweightsadjust(lcdx10$dist,lcdx10$support,lcdx10_pv$support,lcdx10_pv$pv))
-clipw(interpweightsadjust(hy10_5y$dist,hy10_5y$support,hy10_5y_pv$support,hy10_5y_pv$pv))
-clipw(interpweightsadjust(hy10_7y$dist,hy10_7y$support,hy10_7y_pv$support,hy10_7y_pv$pv))
-clipw(interpweightsadjust(lcdx12$dist,lcdx12$support,lcdx12_pv$support,lcdx12_pv$pv))
|