Last data update: 2014.03.03

R: Kernel density estimation for heaped data
dheapingR Documentation

Kernel density estimation for heaped data

Description

Kernel density estimation for heaped data

Usage

dheaping(xheaped, rounds, burnin = 5, samples = 10, setBias = FALSE,
  weights = NULL, bw = "nrd0", boundary = FALSE, unequal = FALSE,
  random = FALSE, adjust = 1, recall = F, recallParams = c(1/3, 1/3))

Arguments

xheaped

heaped values from which to estimate density of x

rounds

rounding values, numeric vector of length >=1

burnin

burn-in sample size

samples

sampling iteration size

setBias

if TRUE a rounding Bias parameter is estimated. For values above 0.5, the respondents are more prone to round down, while for values < 0.5 they are more likely to round up

weights

optional numeric vector of sampling weights

bw

bandwidth selector method, defaults to "nrd0" see density for more options

boundary

TRUE for positive only data (no positive density for negative values)

unequal

if TRUE a probit model is fitted for the rounding probabilities with log(true value) as regressor

random

if TRUE a random effect probit model is fitted for rounding probabilities

adjust

as in density, the user can multiply the bandwidth by a certain factor such that bw=adjust*bw

recall

if TRUE a recall error is introduced to the heaping model

recallParams

recall error model parameters expression(nu) and expression(eta). Default is c(1/3, 1/3)

Value

The function returns a list object with the following objects (besides all input objects):

meanPostDensity

Vector of Mean Posterior Density

gridx

Vector Grid on which density is evaluated

resultDensity

Matrix with Estimated Density for each iteration

resultRR

Matrix with rounding probability threshold values for each iteration (on probit scale)

resultBias

Vector with estimated Bias parameter for each iteration

resultBeta

Vector with estimated Beta parameter for each iteration

resultX

Matrix of true latent values X estimates

Examples

#Simple Rounding   ----------------------------------------------------------
xtrue=rnorm(3000)
xrounded=round(xtrue)
est <- dheaping(xrounded,rounds=1,burnin=20,samples=50)
plot(est,trueX=xtrue)

#####################
#####Heaping
#####################
#Real Data Example  ----------------------------------------------------------
# Student learning hours per week
data(students)
xheaped <- as.numeric(na.omit(students$StudyHrs))
## Not run: est <- dheaping(xheaped,rounds=c(1,2,5,10), boundary=TRUE, unequal=TRUE,burnin=20,samples=50)
plot(est)
summary(est)
## End(Not run)

#Simulate Data   ----------------------------------------------------------
Sim1 <- createSim.Kernelheaping(n=500, distribution="norm",rounds=c(1,10,100),
thresholds=c(-0.5244005, 0.5244005), sd=100)
## Not run: est <- dheaping(Sim1$xheaped,rounds=Sim1$rounds)
plot(est,trueX=Sim1$x)
## End(Not run)

#Biased rounding
Sim2 <- createSim.Kernelheaping(n=500, distribution="gamma",rounds=c(1,2,5,10),
                     thresholds=c(-1.2815516, -0.6744898, 0.3853205),downbias=0.2,
                     shape=4,scale=8,offset=45)
## Not run: est <- dheaping(Sim2$xheaped, rounds=Sim2$rounds, setBias=T, bw="SJ")
plot(est, trueX=Sim2$x)
summary(est)
tracePlots(est)
## End(Not run)

Sim3 <- createSim.Kernelheaping(n=500, distribution="gamma",rounds=c(1,2,5,10),
thresholds=c(1.84, 2.64, 3.05), downbias=0.75, Beta=-0.5, shape=4, scale=8)
## Not run: est <- dheaping(Sim3$xheaped,rounds=Sim3$rounds,boundary=TRUE,unequal=TRUE,setBias=T)
plot(est,trueX=Sim3$x)
## End(Not run)

Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(Kernelheaping)
Loading required package: MASS
Loading required package: ks
Loading required package: KernSmooth
KernSmooth 2.23 loaded
Copyright M. P. Wand 1997-2009
Loading required package: misc3d
Loading required package: mvtnorm
Loading required package: rgl
Loading required package: sparr
Loading required package: spatstat
Loading required package: nlme
Loading required package: rpart

spatstat 1.45-2       (nickname: 'Caretaker Mode') 
For an introduction to spatstat, type 'beginner' 


Attaching package: 'spatstat'

The following object is masked from 'package:MASS':

    area



Welcome to 'sparr': SPAtial Relative Risk (v0.3-8)
T.M. Davies, M.L. Hazelton & J.C. Marshall
-type 'help("sparr")' for details
-type 'citation("sparr")' for how to cite use of this package

> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/Kernelheaping/dheaping.Rd_%03d_medium.png", width=480, height=480)
> ### Name: dheaping
> ### Title: Kernel density estimation for heaped data
> ### Aliases: dheaping
> 
> ### ** Examples
> 
> #Simple Rounding   ----------------------------------------------------------
> xtrue=rnorm(3000)
> xrounded=round(xtrue)
> est <- dheaping(xrounded,rounds=1,burnin=20,samples=50)
[1] "Iteration: 1 of 70"
[1] "Iteration: 2 of 70"
[1] "Iteration: 3 of 70"
[1] "Iteration: 4 of 70"
[1] "Iteration: 5 of 70"
[1] "Iteration: 6 of 70"
[1] "Iteration: 7 of 70"
[1] "Iteration: 8 of 70"
[1] "Iteration: 9 of 70"
[1] "Iteration: 10 of 70"
[1] "Iteration: 11 of 70"
[1] "Iteration: 12 of 70"
[1] "Iteration: 13 of 70"
[1] "Iteration: 14 of 70"
[1] "Iteration: 15 of 70"
[1] "Iteration: 16 of 70"
[1] "Iteration: 17 of 70"
[1] "Iteration: 18 of 70"
[1] "Iteration: 19 of 70"
[1] "Iteration: 20 of 70"
[1] "Iteration: 21 of 70"
[1] "Iteration: 22 of 70"
[1] "Iteration: 23 of 70"
[1] "Iteration: 24 of 70"
[1] "Iteration: 25 of 70"
[1] "Iteration: 26 of 70"
[1] "Iteration: 27 of 70"
[1] "Iteration: 28 of 70"
[1] "Iteration: 29 of 70"
[1] "Iteration: 30 of 70"
[1] "Iteration: 31 of 70"
[1] "Iteration: 32 of 70"
[1] "Iteration: 33 of 70"
[1] "Iteration: 34 of 70"
[1] "Iteration: 35 of 70"
[1] "Iteration: 36 of 70"
[1] "Iteration: 37 of 70"
[1] "Iteration: 38 of 70"
[1] "Iteration: 39 of 70"
[1] "Iteration: 40 of 70"
[1] "Iteration: 41 of 70"
[1] "Iteration: 42 of 70"
[1] "Iteration: 43 of 70"
[1] "Iteration: 44 of 70"
[1] "Iteration: 45 of 70"
[1] "Iteration: 46 of 70"
[1] "Iteration: 47 of 70"
[1] "Iteration: 48 of 70"
[1] "Iteration: 49 of 70"
[1] "Iteration: 50 of 70"
[1] "Iteration: 51 of 70"
[1] "Iteration: 52 of 70"
[1] "Iteration: 53 of 70"
[1] "Iteration: 54 of 70"
[1] "Iteration: 55 of 70"
[1] "Iteration: 56 of 70"
[1] "Iteration: 57 of 70"
[1] "Iteration: 58 of 70"
[1] "Iteration: 59 of 70"
[1] "Iteration: 60 of 70"
[1] "Iteration: 61 of 70"
[1] "Iteration: 62 of 70"
[1] "Iteration: 63 of 70"
[1] "Iteration: 64 of 70"
[1] "Iteration: 65 of 70"
[1] "Iteration: 66 of 70"
[1] "Iteration: 67 of 70"
[1] "Iteration: 68 of 70"
[1] "Iteration: 69 of 70"
[1] "Iteration: 70 of 70"
> plot(est,trueX=xtrue)
> 
> #####################
> #####Heaping
> #####################
> #Real Data Example  ----------------------------------------------------------
> # Student learning hours per week
> data(students)
> xheaped <- as.numeric(na.omit(students$StudyHrs))
> ## Not run: 
> ##D est <- dheaping(xheaped,rounds=c(1,2,5,10), boundary=TRUE, unequal=TRUE,burnin=20,samples=50)
> ##D plot(est)
> ##D summary(est)
> ## End(Not run)
> 
> #Simulate Data   ----------------------------------------------------------
> Sim1 <- createSim.Kernelheaping(n=500, distribution="norm",rounds=c(1,10,100),
+ thresholds=c(-0.5244005, 0.5244005), sd=100)
> ## Not run: 
> ##D est <- dheaping(Sim1$xheaped,rounds=Sim1$rounds)
> ##D plot(est,trueX=Sim1$x)
> ## End(Not run)
> 
> #Biased rounding
> Sim2 <- createSim.Kernelheaping(n=500, distribution="gamma",rounds=c(1,2,5,10),
+                      thresholds=c(-1.2815516, -0.6744898, 0.3853205),downbias=0.2,
+                      shape=4,scale=8,offset=45)
> ## Not run: 
> ##D est <- dheaping(Sim2$xheaped, rounds=Sim2$rounds, setBias=T, bw="SJ")
> ##D plot(est, trueX=Sim2$x)
> ##D summary(est)
> ##D tracePlots(est)
> ## End(Not run)
> 
> Sim3 <- createSim.Kernelheaping(n=500, distribution="gamma",rounds=c(1,2,5,10),
+ thresholds=c(1.84, 2.64, 3.05), downbias=0.75, Beta=-0.5, shape=4, scale=8)
> ## Not run: 
> ##D est <- dheaping(Sim3$xheaped,rounds=Sim3$rounds,boundary=TRUE,unequal=TRUE,setBias=T)
> ##D plot(est,trueX=Sim3$x)
> ## End(Not run)
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>