Datasets

data = iris
datatable(iris)

Dari data tersebut akan dibentuk sebuah scatter plot yang menggunakan fungsi *ggplot* dan *qplot* dengan x =Sepal.Width, y = Petal.Length. Beberapa argument yang dgunakan pada fungsi ggplot adalah

• aes (yang berisiskan variabel yang ingin di deskripsikan dalam hal ini kita dapat mengisikan x =Sepal.Width, y = Petal.Length

• labs yang berisikan judul dari plot

• xlab adalah nama dari variabel

• ylab adalah nama variabel y

Scatterplot

Scatter plot merupakan suatu diagram yang digunakan untuk melihat hubungan dari nilai x dan y berdasarkan sebaran titiknya.Pada data simulasi berikut x adalah nilai siswa dan y adalah persentase dari kehadiran siswa. Berikut diberikan dua cara untuk menampilkan scatterplot.

sc <- ggplot(iris)+
    geom_point(aes(x=Sepal.Width,y = Petal.Length))+
                   labs(title = "Scatter plot")  + 
                   xlab (" Sepal Width")+
                   ylab ("Petal.Length")

# atau dengan 
sc1 <- qplot(x = data$Sepal.Width,y =data$Petal.Length, 
             geom = "point",
             colour = "pink",
             main = "Scater Plot",
             xlab = "Sepal Width",
             ylab = "Petal Lenght")
sc1

Menamppilkan beberapa sebaran data dengan kriteria tertentu dengan menggunakan scatterplot dimana x =Sepal.Width, y = Petal.Length dengan kriteria atau color = Species yang disajikan dapat berbentuk sebagai berikut:

ggplot(iris) +
  geom_point(aes(x =Sepal.Width, y = Petal.Length, color = Species)) + 
  facet_wrap(~Species) + 
  theme(legend.position = "none")

Histogram

data <- iris
head(data)

##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa

ggplot(data) +
    geom_histogram(aes(x = Petal.Length), fill="lightblue",binwidth = 1)+
    ggtitle("Histogram")+
    ylab("Sepal Widht")+
    xlab("Petal Lenght")

Barplot

ggplot(data)+
  geom_bar(aes(x = Species, y = Petal.Length),fill = "lightblue",binwidth = 1,stat = "identity") + 
  ggtitle("Bar Plot")+
    ylab("Sepal Widht")+
    xlab("Petal Lenght")

## Warning: Ignoring unknown parameters: binwidth

ggplot(iris) +
  geom_bar(aes(x =Sepal.Width, y = Petal.Length, color = Species,fill = Species),stat = "identity")+
  facet_wrap(~Species) + 
  theme(legend.position = "none")

Inverse Transform Method

Suatu metode pembangkitan bilangan acak dengan menggunakan inverse dari fungsi cdfnya. tetapi perlu diperhatikan bahwa terdapat kesulitan apabila untuk beberapa sebaran yang tidak dapat di tentukan cdfnya. Namun, metode ini merupakan metode yang fleksibel bahkan dapat digukan pada sebaran diskret.

##  Inverse-Transform method
set.seed(1501211039)
r <- runif(100)
x <- qbeta(r, 3, 2)  
hist(x, prob = TRUE,col ="lightblue") 
sbx <- seq(7, 14, 0.01)  
lines(sbx, dnorm(sbx, mean = 10, sd = 1), col ="blue")  

head(x)

## [1] 0.7266986 0.5845417 0.7052912 0.2634162 0.3444652 0.4695590

set.seed(!501211039)
empirik <- rbeta(100,3,2)
hist(empirik,main="Histogram fungsi Beta",col="lightblue")

Acceptance Rejection Method

Suatu metode ini dapat menjadi suatu alternatif ketika adanya kendala pada integral atau cdf metode Inverse-Transform. Acceptance-rejection method membangkitkan sutu sebaran dari sebaran lainya yang memiliki kriteria sebaran dan panjang data sama dengan data yang ingin kita bangkitkan.

Urutan yang dilakukan:

• Tetapkan peubah acak Y dengan fungsi peluang g yang memiliki suatu konstanta dengan nilai maximum dari \(f(t)/g(t) ≤ c\). Dengan nilai \(f(t) > 0\)

• Bangkitkan u dari sebaran \(Uniform (0,1)\)

• Hitunglah nilai dari \(c * g(y) * u < f(y)\) Apabila terpenuhi maka terima y sebagai x, abila tidak terpenuji mana ulangi langkah 2(a).

Contoh 1 : Membangkitkan suatu fungsi fkp \(6 x (1 - x) ~ beta(2,2)\)

x <- seq(0, 1, 0.01)
y <- 12 * x^2 * (1 - x)
max(y)

## [1] 1.777644

set.seed(1501211039)
n <- 100
j <- k <- 0
y <- numeric(n)
while (k < n) {
    u <- runif(1)
    j <- j + 1
    x <- runif(1)
    if ((12 * x^2 * (1 - x)*dunif(x, 0, 1))> u) {
        k <- k + 1
        y[k] <- x
    }
}
hist(y, prob = T,col = "lightblue")
sbx <- seq(0, 1, 0.01)
lines(sbx, dbeta(sbx, 3, 2),col = "red")

j

## [1] 141

Contoh 2 : Membangkitkan suatu fungsi fkp \(f(x) = 2x^3 +x\) yang dianggap memiiki suatu distribusi tertentu

curve(2*(x^3)+x,0,1)

x <- seq(0, 1, 0.01)
y <- (2*x^3)+x
max(y)

## [1] 3

set.seed(234)
n <- 1000
j <- 0
k <- 0
y <- numeric(n)
while (k < n) {
    u <- runif(1)
    j <- j + 1
    x <- runif(1)
    if (u < (2*(x^3)+x)/3.5 * dunif(x, 0, 1)) {
        k <- k + 1
        y[k] <- x
    }
}
hist(y, prob = T,col = "lightblue")
sbx <- seq(0, 1, 0.01)
lines(sbx, ((2 * (sbx^3))+ x),col ="red")

j

## [1] 3563

Direct Transformation

Metode ini meruupakan membangkitkan sebaran acak berdasarkan transformasi dari sebaran lain. contoh kita ingin membangitkan sebaran khi-kuarat dengan db = 1 hal ini dapat dilakuakn dengan cara membangkitkan suatu sebaran X yang merupakan sebaran normal dengan transformasi \(Y = X^2\).

###  Direct Transformation

set.seed(234)
y <- rnorm(1000)
x <- y^2
hist(x, prob = T,col="lightblue")
sbx <- seq(0, 13, 0.01)
lines(sbx, dchisq(sbx, df = 1), col = "red")

Membangkitkan sebaran T-student db = 4, dengan membangkitkan suatu sebaran \(Khi-kuadrat db =4\) dan \(Normal(0,1)\) dengan transformasi \(Y = z/(sqrt(chi/4)\)

set.seed(454)
z <- rnorm(1000,0,1)
chi <- rchisq(1000,4)
x <- z/(sqrt(chi/4))
hist(x, prob = T,col="lightblue")
sbx <- seq(-5, 5, 0.01)
dt <- dt(sbx, df = 4)
lines(sbx,dt)

Membangkitkan sebaran \(Z ~ Beta(a= 3, b =2)\) , dengan membangkitkan \(X1 ~ gamma(a=3,b=1)\) dan \(X2 ~ gamma(a=2,b=1)\) dengan \(Y = U/(u+v)\)

n <- 1000; a <- 3; b <- 2;
u <- rgamma (n , shape =a , rate =1);
v <- rgamma (n , shape =b , rate =1);
x <- u/( u + v )

hist(x, prob = T,col="lightblue")
sbx <- seq(0, 10, 0.01)
dt <- dbeta(sbx,3,2)
lines(sbx,dt)

set.seed(234)
n<-1000
x1 <- rnorm(n,mean=0,sd=1)
x2<- runif(n,0,1)
x3 <-rgamma(n,2,4)
e <- rnorm(n,mean =0,sd=1)
Y <- 2 +1.2*x1-0.32*x2+3.1*x3 + e

reg <- data.frame(Y,x1,x2,x3)
modelRegresi <- lm(Y ~ ., data = reg)
summary(modelRegresi)

## 
## Call:
## lm(formula = Y ~ ., data = reg)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.2813 -0.6860  0.0271  0.6615  3.9402 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.13768    0.07638  27.987  < 2e-16 ***
## x1           1.20831    0.03271  36.943  < 2e-16 ***
## x2          -0.48967    0.11100  -4.412 1.14e-05 ***
## x3           3.07261    0.09274  33.133  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.017 on 996 degrees of freedom
## Multiple R-squared:  0.7133, Adjusted R-squared:  0.7125 
## F-statistic: 826.1 on 3 and 996 DF,  p-value: < 2.2e-16