Advertisement
simondp

Ellipse BRMS

Nov 17th, 2019
18
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
R 1.36 KB | None
  1. ### Ellipse simulation
  2. ### Simulate X ellipses with y points and attempt to estimate the shape
  3.  
  4. # Load library
  5. library("brms")
  6.  
  7. # Set Seed
  8. set.seed(42)
  9.  
  10. # Define ellipses
  11. nlev <- 8 # Points
  12. nr <- 10 # Participants
  13. width <- 12
  14. height <- 10
  15.  
  16. # Generate n Circles
  17. angs <- seq(pi/2,-pi-2*pi/nlev,-2*pi/nlev)
  18. capture <- matrix(NA,nlev,3)
  19. capture[,1] <- 1
  20. # Add noise per meassure
  21. indvar <- rnorm(1)
  22. random <- rnorm(nlev,indvar,sd=0.1)
  23. capture[,2] <- cos(angs)*(width+random)
  24. capture[,3] <- sin(angs)*(height+random)
  25. captureTemp <- matrix(NA,nlev,3)
  26. for(i in 2:nr){
  27.   indvar <- rnorm(1)
  28.   random <- rnorm(nlev,indvar,sd=0.1)
  29.   captureTemp[,1] <- i
  30.   captureTemp[,2] <- cos(angs)*(width+random)
  31.   captureTemp[,3] <- sin(angs)*(height+random)
  32.   capture <- rbind(capture,captureTemp)}
  33. # Make data frame
  34. dat <- data.frame(capture)
  35. names(dat) <- c("part", "xpos", "ypos")
  36.  
  37. # Plot
  38. limit <- max(capture)
  39. plot(dat$xpos,dat$ypos, xlim = c(-limit,limit), ylim = c(-limit,limit))
  40.  
  41. # Brms (wrong method, surely)
  42. m1 <- brm(data = dat,
  43.                    family = gaussian,
  44.                    mvbind(xpos, ypos) ~ 1 +(1|part), chains = 2, cores = 2, iter=500)
  45. m1
  46. # Attempt at better model
  47. m2 <- brm(data = dat,
  48.              family = gaussian,
  49.              mvbind(xpos*a*(cos(t)),
  50.                     ypos*b*(sin(t))) ~ 1 +(1|part), chains = 2, cores = 2, iter=500)
Advertisement
RAW Paste Data Copied
Advertisement