# Plot logarithm(theta) vs theta

I have created the below function for Em two binomial mixture and i want to plot the log(theta) versus theta of each iteration but i don't know how to do it?Any help i will appreciate it.

``````
EMflipcoin = function(p,n,theta1,theta2,x) {

tol=10^(-10)
maxiter=100
diff=1
iter=0
while(diff>tol & iter<maxiter){
theta1 = theta1
theta2 = theta2
p=p
## E-step:
a   = (theta1^x)*(1-theta1)^(n-x);a
b   = (theta2^x)*(1-theta2)^(n-x);b
pa  = a/(a+b);pa
pb  = b/(a+b);pb
eah = pa*x;eah
eat = pa*(n-x);eat
ebh = pb*x;ebh
ebt = pb*(n-x);ebt
c   = dbinom(x,n,p);c
c2  = dbinom(x,n,1-p);c2
pc  = c/(c+c2);pc

## M-step:
theta1new = sum(eah)/(sum(eah)+sum(eat));theta1new
theta2new = sum(ebh)/(sum(ebh)+sum(ebt));theta2new
pnew      = mean(pc);pnew
## calculate diff for tolerance in while condition
diff = (pnew-p)^2+(theta1new-theta1)^2+(theta2new-theta2)^2
iter =iter+1;

theta1 = theta1new
theta2 = theta2new
p =pnew
cat("Iter", iter,
", theta1=",theta1new,
", theta2=",theta2new,
", p     =",pnew,
", diff=", diff, "\n")}

}
x=c(5,9,8,4,7)
EMflipcoin(p=0.5,n=10,theta1=0.7,theta2=0.4,x=x)

``````

Which one of the thetas you are interested in to plot? Anyway, you need to store the initial value in some variables and add the new values for each iteration. Then, add the plot instructions at the end of your function. Here is a proposal using theta1:

EMflipcoin = function(p,n,theta1,theta2,x) {
tol=10^(-10)
maxiter=100
diff=1
iter=0

``````value_theta = theta1 # save original values
value_log_theta = log(theta1)

while(diff>tol & iter<maxiter){
theta1 = theta1
theta2 = theta2
p=p

## E-step:
a   = (theta1^x)*(1-theta1)^(n-x);a
b   = (theta2^x)*(1-theta2)^(n-x);b
pa  = a/(a+b);pa
pb  = b/(a+b);pb
eah = pa*x;eah
eat = pa*(n-x);eat
ebh = pb*x;ebh
ebt = pb*(n-x);ebt
c   = dbinom(x,n,p);c
c2  = dbinom(x,n,1-p);c2
pc  = c/(c+c2);pc

## M-step:
theta1new = sum(eah)/(sum(eah)+sum(eat));theta1new
theta2new = sum(ebh)/(sum(ebh)+sum(ebt));theta2new
pnew      = mean(pc);pnew
## calculate diff for tolerance in while condition
diff = (pnew-p)^2+(theta1new-theta1)^2+(theta2new-theta2)^2
iter =iter+1;

theta1 = theta1new
theta2 = theta2new
p =pnew
cat("Iter", iter,
", theta1=",theta1new,
", theta2=",theta2new,
", p     =",pnew,
", diff=", diff, "\n")

value_theta = c(value_theta,theta1)
value_log_theta = c(value_log_theta, log(theta1) )
}

# plot theta vs log(theta)
df <- data.frame(value_log_theta, value_theta) #create a data frame with your variables
plot <- ggplot(df,aes(x = value_log_theta, y = value_theta)) +
geom_point() +
geom_line() +
labs(x = "log(theta)", y = "theta") +
theme_bw()
print(plot)
``````

}
x=c(5,9,8,4,7)
EMflipcoin(p=0.5,n=10,theta1=0.7,theta2=0.4,x=x)

This topic was automatically closed 21 days after the last reply. New replies are no longer allowed.