[R] what meaning missing value True /False needed
Sarah Goslee
sarah.goslee at gmail.com
Fri Jul 31 15:21:50 CEST 2009
You are far more likely to get helpful volunteers to answer your
questions if you
distill your problem down to the smallest point. That's a lot of code to wade
through. But since I was waiting for something else to run, I took a stab at it.
It also would have been helpful if you'd given the *complete* error message,
since it tells you where exactly the problem is.
In this case, the problem is in LapseRate(), where PriceFactor is compared to
WindowLeft (as the actual error message says). I added two lines using cat()
to print those two variables each time LapseRate() is called, and here's what
I saw:
Maximization Problem.
PriceFactor: 0.9352139
WindowLeft: 1
PriceFactor: 1.627354
WindowLeft: 1
PriceFactor: NA
WindowLeft: 1
Error in if (PriceFactor < WindowLeft) r <- CentralLapseRate * (1 +
ElasticityLeft/WindowLeft * :
missing value where TRUE/FALSE needed
So indeed, LapseRate() is comparing a NA value to something else,
resulting in NA,
just as the error message stated.
Figuring out why PriceFactor is NA after two iterations is left as an
exercise for the
reader.
Sarah
On Fri, Jul 31, 2009 at 8:21 AM, Inchallah Yarab<inchallahyarab at yahoo.fr> wrote:
> This is my code i don't understand the error message:
>
> library(rgenoud)
> rm(list=ls())
> set.seed(666)
> #########################################################
> # As a first step, it is assumed that all input parameters are independent of ageing :
> #########################################################
> InputDim <-20
> # Max number of ageings in the inputs
> CPIRate <- rep(0.02 , InputDim )
> # CPI inflation rate
> LossInflaRate <- rep(0.02 , InputDim)
> # inflation rate of losses
> RiskFreeRate <- rep(0.04, InputDim )
> # vector of risk free rate
> RiskRate <- rep(0.06, InputDim )
> # vector of risk free rate
> ################ lapse param begin :
> LapseType <- rep("linear" , InputDim )
> # type of the function used to calculate lapse rate= "power", "logistic", "linear"
> CentralLapseRate <- rep(0.13 , InputDim )
> # lapse rate for price factor = 1
> ElasticityLeft <- rep( 10, InputDim )
> # left elasticity for price factor = 1
> ElasticityRight <- rep(10, InputDim )
> # right elasticity for price factor = 1
> WindowLeft <- rep( 1 , InputDim )
> # friction window price factor left bound
> WindowRight <- rep( 1, InputDim )
> # friction window price factor right bound
> MinLapseRate <- rep( 0 , InputDim )
> # min lapse rate
> MaxLapseRate <- rep( 1 , InputDim )
> # min lapse rate
> Alpha <- 0.5
> #acquisition operator
> ############### lapse param end
> Price <- rep(1 , InputDim )
> # Axa price for central market position and price factor = 1
> IndivMeanLossAmount <- rep(0.97 , InputDim )
> # mean individual loss amount
> Commission <- rep(0 , InputDim )
> # commision rate
> IndividualExp <- rep(0 , InputDim )
> # expense amount per risk
> #######################################
> # Market Expert Model parameter begin :
> #######################################
> LRUpBound <- 1.01349495023810
> # Upper bound of the loss ratio
> LRLowBound <- 0.743806254238095
> # Lower bound of the loss ratio
> TransitCoefMean <- 1.03001118178481
> # Mean of the coefficient of transition
> TransitCoefStdDev <- 0.0187801387384363
> # Standard deviation of the coefficient of transition
> #######################################
> # Optimization algorithm & Monte-Carlo parameters :
> #######################################
> LTMeanComputYearNb <- 1000
> # Number of year used to simulate LR
> NbOfSimScenar <- 1000
> # Number Monte-Carlo scenarios
> GenAlgoIterNb <- 100
> # number of iteration for genetic algorithm
> #######################################
> # Problem size parameters :
> #######################################
> FirstYear <-2009
> # first year of renewal
> LastYear <- 2018
> # last year of projected renewal
> ######################################
> # Function giving the minimum x (price factor = 1+x ) :
> ######################################
> Xmin <- function(LapseType, ElasticityLeft, WindowLeft, MinLapseRate ,CentralLapseRate)
> {
> InputDim <- length(LapseType)
> x=numeric(InputDim)
> for (i in 1:InputDim)
> {
> if (LapseType[i]=="power" )
> x[i] <- WindowLeft [i]*(MinLapseRate[i] /CentralLapseRate[i])^(1/ElasticityLeft[i])-1
> else if (LapseType[i]=="linear")
> x [i]<- (WindowLeft[i]/ElasticityLeft[i])*(ElasticityLeft[i]-1+MinLapseRate[i] /CentralLapseRate[i])-1
> }
> x
> }
>
> ######################################
> # Function giving the maximum x (price factor = 1+x ) :
> ######################################
> Xmax <- function(LapseType, ElasticityLeft, WindowLeft, MaxLapseRate ,CentralLapseRate)
> {
> InputDim <- length(LapseType)
> x=numeric(InputDim)
> for (i in 1:InputDim)
> {
> if (LapseType[i]=="power" )
> x[i] <- WindowLeft[i] *(MaxLapseRate[i]/CentralLapseRate[i])^(1/ElasticityLeft[i])-1
> else if (LapseType[i]=="linear")
> x[i] <- (WindowLeft[i]/ElasticityLeft[i])*(ElasticityLeft[i]-1+MaxLapseRate[i]/CentralLapseRate[i])-1
> }
> x
> }
>
> ##########################################################
> # Function simulating a scenario of projected market loss ratio :
> ##########################################################
> LRScenar <- function(PriorLR,Phase,LRUpBound ,LRLowBound,TransitCoefMean,TransitCoefStdDev,LTMeanComputYearNb)
> {
> initial <- PriorLR
> LRScenar <- NULL
> phase2 <- Phase
> SigmaLog <- (log(1+(TransitCoefStdDev /TransitCoefMean)^2))^0.5
> MuLog <- log(TransitCoefMean)-SigmaLog^2/2
> for (i in 1:LTMeanComputYearNb)
> {
> kk <- rlnorm(1, meanlog = MuLog, sdlog = SigmaLog)
> if (phase2!=1) kk <- (1/kk)
> x_1 <- initial*kk
> initial <- x_1
> if (((x_1 > LRUpBound )&(phase2==1))|((x_1<LRLowBound)&(phase2==-1))) phase2 <- (-phase2)
> LRScenar <- c(LRScenar,x_1)
> }
> LRScenar
> }
> ################################################
> # Function giving the lapse rate :
> ################################################
> LapseRate <- function(x, LapseType,CentralLapseRate,
> ElasticityLeft, ElasticityRight,WindowLeft, WindowRight,
> MinLapseRate,MaxLapseRate )
> {
> PriceFactor <- 1+x
> if (LapseType=="power")
> {
> if(PriceFactor<WindowLeft)
> r <- CentralLapseRate*(PriceFactor/WindowLeft)^ElasticityLeft
> else if((PriceFactor<=WindowRight)&(PriceFactor>=WindowLeft ))
> r <- CentralLapseRate
> else
> r <-CentralLapseRate*(PriceFactor/WindowRight)^ElasticityRight
> }
> else if (LapseType=="linear")
> {
> if(PriceFactor<WindowLeft)
> r <- CentralLapseRate*(1+ElasticityLeft/WindowLeft*(PriceFactor-WindowLeft))
> else if((PriceFactor<=WindowRight)&(PriceFactor>=WindowLeft ))
> r <- CentralLapseRate
> else
> r <- CentralLapseRate*(1+ElasticityRight/WindowRight*(PriceFactor-WindowRight))
> }
> min( max( r, MinLapseRate ) , MaxLapseRate)
> }
>
> ######################################################
> # Function calculating a matrix of market price (excluding inflation) scenarios :
> ######################################################
> MarketPriceMatrix <- function(NbOfFutureSimYear ,NbOfSimScenar,Phase,LRUpBound ,LRLowBound,TransitCoefMean,
> TransitCoefStdDev,PriorMarketPricePosition,MarketPricePhase,MarketPriceLTMean)
> {
> LRPhase <- -MarketPricePhase
> MarketPriceMatrix <- matrix(numeric(NbOfFutureSimYear*NbOfSimScenar), nrow=NbOfSimScenar,ncol=NbOfFutureSimYear)
> PriorLR <- 1 / (MarketPriceLTMean * PriorMarketPricePosition)
> for (i in 1:NbOfSimScenar)
> MarketPriceMatrix[i,] <- 1/ LRScenar (PriorLR,LRPhase,LRUpBound ,LRLowBound,TransitCoefMean,
> TransitCoefStdDev,NbOfFutureSimYear)/MarketPriceLTMean
> MarketPriceMatrix
> }
> #############################################
> # Function calculating a scenario of results :
> #############################################
> ResultScenar <- function(FirstYear,OptimYear,Price,CumulRetentionCoef,Commission,IndividualExp,NbOfFutureSimYear,
> IndivMeanLossAmount,CPIRate,LossInflaRate,Alpha)
> {
> CumulLossInflaRate <- cumprod(1+LossInflaRate)
> CumulCPIRate <- cumprod(1+CPIRate)
> y <- numeric(NbOfFutureSimYear)
> delta <- OptimYear - FirstYear
> if (OptimYear==FirstYear)
> {
> y[1] <CumulRetentionCoef[1]*( (Price[1]*(1-Commission[delta +1])-IndivMeanLossAmount[delta +1])*CumulLossInflaRate[delta +1]
> - IndividualExp[delta +1]*CumulCPIRate[delta +1])
> }
> for (i in 2: NbOfFutureSimYear)
> {
> k <- delta + i
> y[i] <- (1-Alpha)*CumulRetentionCoef[i-1]*( (Price[i-1]*(1-Commission[k-1])-IndivMeanLossAmount[k-1])*CumulLossInflaRate[k-1]
> - IndividualExp[k]*CumulCPIRate[k-1] )
> + Alpha*CumulRetentionCoef[i]*( (Price[i]*(1-Commission[k])-IndivMeanLossAmount[k])*CumulLossInflaRate[k]
> - IndividualExp[k]*CumulCPIRate[k] )
> }
> if (OptimYear!=FirstYear)
> {
> for (i in 1: (NbOfFutureSimYear+1))
> {
> k <- delta + i
> y[i] <- (1-Alpha)*CumulRetentionCoef[i]*( (Price[i]*(1-Commission[k-1])-IndivMeanLossAmount[k-1])*CumulLossInflaRate[k-1]
> - IndividualExp[k]*CumulCPIRate[k-1] )
> + Alpha*CumulRetentionCoef[i+1]*( (Price[i+1]*(1-Commission[k])-IndivMeanLossAmount[k])*CumulLossInflaRate[k]
> - IndividualExp[k]*CumulCPIRate[k] )
> }
> }
> y
> }
> ##############################################
> ### Function calculating a scenario of results :
> ###############################################
> ResultActuSum <- function(FirstYear,OptimYear,x, NbOfFutureSimYear,
> CentralLapseRate, ElasticityLeft,WindowLeft, WindowRight, ElasticityRight,
> MarketPriceMean,Price,Commission,IndividualExp,IndivMeanLossAmount,
> CPIRate,LossInflaRate,Beta,Alpha)
> {
> R <- numeric(NbOfFutureSimYear+1)
> Price2 <- numeric(NbOfFutureSimYear+1)
> CumulRetentionCoef <- numeric(NbOfFutureSimYear+1)
> Result <- numeric(NbOfFutureSimYear+1)
> delta <- OptimYear - FirstYear
> for (i in 1:(NbOfFutureSimYear+1))
> {
> k <- delta + i
> R[i] <- LapseRate(x[i], LapseType[k],CentralLapseRate[k],
> ElasticityLeft[k], ElasticityRight[k],WindowLeft[k], WindowRight[k],
> MinLapseRate[k],MaxLapseRate[k] )
> Price2[i] <- Price[k]*(1+x[i])*MarketPriceMean[i]
> }
> CumulRetentionCoef <- cumprod(1-R)
> Result <- ResultScenar(FirstYear,OptimYear,Price2,CumulRetentionCoef,Commission,IndividualExp,NbOfFutureSimYear, IndivMeanLossAmount,
> CPIRate,LossInflaRate,Alpha)
> sum(Beta[1:NbOfFutureSimYear]*Result)
> }
>
> #########################################################
> # Function calculating the optimal poistioning for a given problem year & filtration :
> #########################################################
> YearFiltrationPbOptimPosition <- function( PriorMarketPricePosition,MarketPricePhase,LRUpBound,LRLowBound,TransitCoefMean,TransitCoefStdDev,
> LTMeanComputYearNb ,
> FirstYear,LastYear,OptimYear,
> CentralLapseRate, ElasticityLeft,WindowLeft, WindowRight, ElasticityRight,
> Price,Commission,IndividualExp,IndivMeanLossAmount,
> CPIRate,LossInflaRate,RiskFreeRate,RiskRate,Alpha)
> {
> Beta <- 1/cumprod(1+RiskFreeRate+RiskRate)
> CumulLossInflaRate <- cumprod(1+LossInflaRate)
> CumulCPIRate <- cumprod(1+CPIRate)
> xmin <- Xmin(LapseType, ElasticityLeft, WindowLeft, MinLapseRate,CentralLapseRate)
> xmax <- Xmax(LapseType, ElasticityLeft, WindowLeft, MaxLapseRate,CentralLapseRate)
> kmin <- OptimYear - FirstYear + 1
> kmax <- LastYear - FirstYear + 1
> NbOfFutureSimYear <- LastYear - OptimYear + 1
> PriorLR <- (LRUpBound*LRLowBound)^0.5
> Phase <- 1
> LRScenar <- LRScenar(PriorLR,Phase,LRUpBound ,LRLowBound,TransitCoefMean,TransitCoefStdDev,LTMeanComputYearNb)
> MarketPriceLTMean <- mean( 1/LRScenar)
> MarketPriceMatrix <- MarketPriceMatrix(NbOfFutureSimYear+1,NbOfSimScenar,Phase,LRUpBound ,LRLowBound,
> TransitCoefMean,TransitCoefStdDev,PriorMarketPricePosition,MarketPricePhase,MarketPriceLTMean)
> MarketPriceMean <- apply(MarketPriceMatrix,2,mean)
> Evaluate <- function(x)
> {
> ResultActuSum(FirstYear,OptimYear,x,NbOfFutureSimYear,
> CentralLapseRate, ElasticityLeft,WindowLeft, WindowRight, ElasticityRight,
> MarketPriceMean,Price,Commission,IndividualExp,IndivMeanLossAmount,
> CPIRate,LossInflaRate,Beta,Alpha)
> }
> output <- genoud(Evaluate , nvars=NbOfFutureSimYear,Domains = cbind(xmin[kmin : kmax],xmax[kmin : kmax]),
> pop.size=5000,max=TRUE)
> Solution <- rev(output$par)
> MarketPriceMean <-rev(MarketPriceMean)
> Solution
> }
> ########################################################
> # Problem parameters for a given year N & the coressponding filtration in year N-1 :
> ########################################################
>
> OptimYear <- 2017
> # year of optimisation
> PriorMarketPricePosition <- 1
> #market price position in OptimYear-1
> MarketPricePhase <- 1
> #market price (excluding inflation) phase in OptimYear
> #1: increasing, -1: decreasing
> ########################################################
> Xoptim <- YearFiltrationPbOptimPosition( PriorMarketPricePosition,MarketPricePhase,LRUpBound,LRLowBound,TransitCoefMean,TransitCoefStdDev,
> LTMeanComputYearNb ,
> FirstYear,LastYear,OptimYear,
> CentralLapseRate, ElasticityLeft,WindowLeft, WindowRight, ElasticityRight,
> Price,Commission,IndividualExp,IndivMeanLossAmount,
> CPIRate,LossInflaRate,RiskFreeRate,RiskRate,Alpha)
>
> Xoptim
>
>
>
> Best Regards
> thank you for your help
>
>
>
> [[alternative HTML version deleted]]
>
>
--
Sarah Goslee
http://www.functionaldiversity.org
More information about the R-help
mailing list