Generating multiple samples from distributions in R
Suppose I wish to generate m samples of n variables of the random variable X and n variables of Y in R.
How would I write such a code? I experimented with using a double forloop with indexing each sample but was unsuccessful doing so.
For ease assume X follows a normal (0,1) distribution and Y a uniform (0,1) distribution.
(Note: I am not particular about the distribution, rather how the code is written)
See also questions close to this topic

Problems with if() condition in Stan/RStan when modelling values from binomial random variable
I'm trying to use Stan and R to fit a model that, uhh, models the observed realisations y_i = 16, 9, 10, 13, 19, 20, 18, 17, 35, 55, which are from a binomial distributed random variable, say, Y_i, with parameters m_i (number of trials) and p_i (success probability in each trial).
yi = c(16, 9, 10, 13, 19, 20, 18, 17, 35, 55)
For the purposes of this experiment, I'm going to assume that all of the m_i are fixed and given by m_i = 74, 99, 58, 70, 122, 77, 104, 129, 308, 119.
mi = c(74, 99, 58, 70, 122, 77, 104, 129, 308, 119)
I'm going to use Jeffrey's prior: \alpha=0.5 and \beta=0.5.
alpha = 0.5, beta = 0.5
I'm trying to
My attempt at 2. is this section of code:
real k; real mx = 0; real mn = 0; if (p > mx) mx = p; if (mn > p) { mn = p; } k = mx  mn;
My Stan code is as follows:
```{stan output.var="BinModBeta"} data { int <lower = 1> mi[10]; int <lower = 0> yi[10]; real <lower = 0> alpha; real <lower = 0> beta; } parameters { real <lower = 0, upper = 1> p[10]; } transformed parameters { real k; real mx = 0; real mn = 0; if (p > mx) mx = p; if (mn > p) { mn = p; } k = mx  mn; } model { yi ~ binomial(mi, p); p ~ beta(alpha, beta); } ```
My R code is as follows:
```{r} library(rstan) ``` ```{r} data.in < list(mi = c(74, 99, 58, 70, 122, 77, 104, 129, 308, 119), yi = c(16, 9, 10, 13, 19, 20, 18, 17, 35, 55), alpha = 0.5, beta = 0.5) model.fit1 < sampling(BinModBeta, data=data.in) ``` ```{r} print(model.fit1, pars = c("p"), probs=c(0.1,0.5,0.9), digits = 5) ```
Now, I just started learning Stan, so I'm honestly not sure if this is correct at all. However, it seems like this code works for my first aim (at least, whatever I've coded seems to work ...). But my trouble starts when attempting to code my second aim.
When I attempt to compile the Stan code above, I get the following error:
Now, based on this error message, it seems like my issue is arising from the fact that p is a vector of 10 real values, rather than a single real. However, due to my inexperience with Stan, I'm completely unsure of how to get around this problem.
I would greatly appreciate it if people could please take the time to help me fix this.

combining JSON pages with RJSONIO
I am having a problem with combining multiple pages with
RJSONIO
.Base version: If I try with one page, it works:
test11<RJSONIO::fromJSON("http://zoeken.kvk.nl/Address.ashx?site=handelsregister&q=01000040") df<ldply(test11, data.frame)
However, when I try to automate the process for multiple pages, then the list stays empty. I used the following code:
baseurl < "http://zoeken.kvk.nl/Address.ashx?site=handelsregister" pages < list() for(i in 01000039:01000066){ mydata < RJSONIO::fromJSON(paste0(baseurl, "&q=", i)) message("Retrieving page ", i) pages[[i+1]] < mydata$resultatenHR } #combine all into one filings < rbind_pages(pages) #check output nrow(filings)
Nothing happens, and the 'pages' list is empty, and the mydata is also empty. I tried using another json package, but with the other packages (jsonlite and rjason) I get an error message:
"Error in open.connection(con, "rb") : Timeout was reached: Connection timed out after 10000 milliseconds"
Thanks in advance for any reply

mirt package  standard error
I am using
mirt
package to make a2PL
model. I was trying to get thestandard error
for each estimated parameter.After using
mirt(mydata, 1, itemtype = '2PL',SE=TRUE)
, i am getting the following results which I am confused which info shows the standard error.FYI,
mydata
is sparse and binary. Also, I already read this page, which I did not get my answer.Results sample for the item 59
$`59` a1 d g u par 36.859 12.266 0 1 CI_2.5 905.576 445.648 NA NA CI_97.5 979.295 421.117 NA NA
Which information shows the standard error?
a1
,d
,g
,oru
.Appreciate!

For loop, running a 1000 times the program, average number of rolls required to get snake eyes in JAVA
So I posted my question earlier and thought I was all good, but I messed up and realized I had understood the question completely wrong. I do not need to calculate the average of snake eyes over 1000 dice rolls, but the average of number of rolls to get a snake eyes, over a 1000 play. I am a little lost in how to accomplish that. I tried this:
public class RollDiceforloop { public static void main(String[] args) { int die1, die2, snakeye, rolls, game; snakeye = 0; die1 = 0; die2 = 0; rolls = 0; for (game = 0; game < 1000; game++) { die1 = (int)(Math.random()*6)+1; die2 = (int)(Math.random()*6)+1; if (die1 != 1 && die2 != 1); { rolls +=1; } if (die1 == 1 && die2 == 1) { snakeye +=1; rolls +=1; } } float average = snakeye / (float) rolls; TextIO.putln(""+snakeye+" Snake Eyes over "+game+" games, with an average of "+average+" rolls required to get a Snake Eye."); } }
But I am not getting the correct result. I am a bit lost on how to accomplish this. Help please?

Why do I get a cast exception in the "for each" and not in the "for"?
I would expect the two fors in the following code to behave the same way but they don't:
ArgumentCaptor<Appender> argumentCaptor = ArgumentCaptor.forClass(Appender.class); Mockito.verify(mockedAppender, times(3)).doAppend(argumentCaptor.capture()); for(int i = 0; i< argumentCaptor.getAllValues().size(); ++i) { System.out.println(((LoggingEvent) argumentCaptor.getAllValues().get(i)).getMessage()); } for(Appender appender : argumentCaptor.getAllValues()) { System.out.println(((LoggingEvent) appender).getMessage()); }
The first correctly prints the messages while the second gives me this exception:
java.lang.ClassCastException: ch.qos.logback.classic.spi.LoggingEvent cannot be cast to ch.qos.logback.core.Appender

For loop with special condition
public int minCompletionTime() { // finds minCompletionTime int time = 1; for(Job j : jobs) { if (j.getStartTime() == 1) { return 1; } } for (Job j : jobs) { //Calculate the minimum completion time if (j.getStartTime() + j.time > time) { time = j.getStartTime() + j.time; } } return time; }
Can someone explain the loop condition? it looks new to for me. Thank you

Multivariate Normal Vectors with Associated Cost (Classification)
I'm working on a classification problem with two populations, and I need to generate new samples to classify. I'm using mvrnorm() to generate new vectors, but I'm wondering if there is any way to get an associated cost with the generated vectors? I know R has a lot of tricks, but I'm not finding anything for this. Any help would be appreciated.

What is equivalent check for Standard deviation in normal distribution in C++ from Python?
I have Python script for Gaussian Normal distribution:
import numpy as np x_mu = 25 x_sigma = 5 size = 1000 x_distribution = np.random.normal(x_mu, x_sigma, size) #i am looking for help ONLY FOR this line below test_distribution = np.std(x_distribution) print (test_distribution)
So, for larger size, test_distribution value will be closer to x_sigma value.
I searched in Google, but could not find solution for C++ FOR the "test_distribution", not entire code. Please, if you know any libs in C++ or ideas, leave in comments or answer. Thanks

1D histogram and Gaussian fit
Let's consider we have a 1D histogram, and we can fit Gaussian on it. And the fit is reasonably good. But there are always few counts outside the fit.
Now is there any way to subtract the fitted curve from the actual histogram, so that we can get to know the goodness of the fitting?

Monte Carlo standard deviation python
Hello I tried to solve a problem in python with the Monte Carlo method. The problem is the following: A company services notebooks. A review of its records shows that the time taken for a service call is normally distributed with a mean of 60 minutes and standard deviation of 20 minutes. a. What proportion of service calls take less than one hour? b. What proportion of service calls take more than 50 minutes? c. What proportion of service calls take more than 80 minutes? d. If a service call has taken more than an hour, what is the probability that it will take more than 70 minutes? d. In a random sample of four calls, what is the probability that the average length of calls is less than one 50 minutes?
My solution for this:
import math import random calls = 4 sa = 0 sb = 0 sc = 0 sd = 0 sd2 = 0 se = 0 sumMin = 0 for i in range (1,calls): r = 60 + 20.* random.randint(1, calls) #print(r) if r<60: sa = sa + 1 break if r > 50: sb = sb + 1 break if r > 80: sc = sc + 1 break if r > 60: sd2 = sd2 + 1 break if r > 60: sd = sd + 1 break break sumMin = sumMin + r meanCallsl = sumMin / calls if meanCallsl < 50: se = se + 1
and for a,b,c,d the following conditions
pa = 100 * sa/ calls print (pa) pb = 100 * sb / calls print(sb) pc = 100 * sd / calls print(pc) pd = 100*sd/sd2 print (pd) pe = 100 * se/calls print(pe)
Though I am not sure of the outcome, could you please advise?
Thank you

R programming (Monte Carlo Simulation)
Hi I would like to ask how to sample out 50 instances of iris data (which contain 150 instances) by using Monte Carlo Simulation ? Any idea? Many thanks

How to plot a point in some color based on the result of a comparaison, not using a loop?
I'm using random points to determine the area below a curve (MonteCarlo):
 X: 1xn vector of x values for the function
 Y: 1xn vector of y = f(x)
 RP: mxn matrix of m random y for each x
I would like to split RP into RPA and RPB depending on it being above or below the curve. The idea is then to plot RPA and RPB against X, in different colors. This code doesn't work because RPA and RPB number of columns is not the same than X:
clf f = @(x) sin(x/10) + cos(x/60); % Function xMin = 1; xMax = 100; % x interval X = [xMin:xMax]; Y = f(X); plot(X,Y), hold on % Plot function yMin = min(Y); yMax = max(Y); % Axes limits set(gca, 'xlim', [xMin, xMax], 'ylim', [yMin, yMax]) m = 20; % Random points per x value RP = rand(m, columns(X)) .* (yMaxyMin) .+ yMin; % Split points (doesn't work) RPA = RP(RP>Y); RPB = RP(RP<=Y); br = size(RPB) / size(RP) % Ratio of points below a = (xMax  xMin) * (yMax  yMin) * br % Area below % Plot points plot(X, RPA, 'r.') % Above plot(X, RPB, 'g.') % Below
Is there a possibility to build RPA and RPB so that they are the same size than RP, with the excluded points y being NaN or something similar, which can be counted and plotted?