How to generate random integers between 1 and 4 that have a specific mean? I need to generate 100 random integers in R, where each integer is between 1 and 4 (hence 1,2,3,4) and the mean is equal to a specific value.
If I draw random uniform numbers between 1 and 5 and get floor, I have a mean of 2.5.
x = floor(runif(100,min=1, max=5)) 

I need to fix the mean to 1.9 or 2.93 for example.
I guess I can generate random integers that add to 100 * mean but I don't know how to restrict to random integers between 1 and 4.
 A: I agree with X'ian that the problem is under-specified.  However, there is an elegant, scalable, efficient, effective, and versatile solution worth considering.
Because the product of the sample mean and sample size equals the sample sum, the problem concerns generating a random sample of $n$ values in the set $\{1,2,\ldots, k\}$ that sum to $s$ (assuming $n \le s \le kn,$ of course).
To explain the proposed solution and, I hope, justify the claim of elegance, I offer a graphical interpretation of this sampling scheme.  Lay out a grid of $k$ rows and $n$ columns.  Select every cell in the first row.  Randomly (and uniformly) select $s-n$ of the remaining cells in rows $2$ through $k.$ The value of observation $i$ in the sample is the number of cells selected in column $i:$

This $4\times 100$ grid is represented by black dots at the unselected cells and colored patches at the selected cells.  It was generated to produce a mean value of $2,$ so $s=200.$  Thus, $200-100=100$ cells were randomly selected among the top $k-1=3$ rows.  The colors represent the numbers of selected cells in each column. There are $28$ ones, $47$ twos, $22$ threes, and $3$ fours.  The ordered sample corresponds to the sequence of colors from column $1$ through column $n=100.$
To demonstrate scalability and efficiency, here is an R command to generate a sample according to this scheme.  The question concerns the case $k=4, n=100$ and $s$ is $n$ times the desired average of the sample:
tabulate(sample.int((k-1)*n, s-n) %% n + 1, n) + 1

Because sample.int requires $O(s-n)$ time and $O((k-1)n)$ space, and tabulate requires $O(n)$ time and space, this algorithm requires $O(\max(s-n,n))$ time and $O(kn)$ space: that's scalable.  With $k=4$ and $n=100$ my workstation takes only 12 microseconds to perform this calculation: that's efficient.
(Here's a brief explanation of the code.  Note that integers $x$ in $\{1,2,\ldots, (k-1)n\}$ can be expressed uniquely as $x = nj + i$ where $j \in \{0,1,\ldots, k-2\}$ and $i\in\{1,2,\ldots, n\}.$  The code takes a sample of such $x,$ converts them  to their $(i,j)$ grid coordinates, counts how many times each $i$ appears (which will range from $0$ through $k-1$) and adds $1$ to each count.)
Why can this be considered effective?  One reason is that the distributional properties of this sampling scheme are straightforward to work out:


*

*It is exchangeable: all permutations of any sample are equally likely.

*The chance that the value $x \in\{1,2,\ldots, k\}$ appears at position $i,$ which I will write as $\pi_i(x),$ is obtained through a basic hypergeometric counting argument as $$\pi_i(x) = \frac{\binom{k-1}{x-1}\binom{(n-1)(k-1)}{s-n-x+1}}{\binom{n(k-1)}{ s-n}}.$$  For example, with $k=4,$ $n=100,$ and a mean of $2.0$ (so that $s=200$) the chances are $\pi = (0.2948, 0.4467, 0.2222, 0.03630),$ closely agreeing with the frequencies in the foregoing sample.  Here are graphs of $\pi_1(1), \pi_1(2), \pi_1(3),$ and $\pi_1(4)$ as a function of the sum:


*The chance that the value $x$ appears at position $i$ while the value $y$ appears at position $j$ is similarly found as $$\pi_{ij}(x,y) = \frac{\binom{k-1}{x-1}\binom{k-1}{y-1}\binom{(n-1)(k-1)}{s-n-x-y+2}}{\binom{n(k-1)}{ s-n}}.$$
These probabilities $\pi_i$ and $\pi_{ij}$ enable one to apply the Horvitz-Thompson estimator to this probability sampling design as well as to compute the first two moments of the distributions of various statistics.
Finally, this solution is versatile insofar as it permits simple, readily-analyzable variations to control the sampling distribution.  For instance, you could select cells on the grid with specified but unequal probabilities in each row, or with an urn-like model to modify the probabilities as sampling proceeds, thereby controlling the frequencies of the column counts.
A: You can use sample() and select specific probabilities for each integer. If you sum the product of the probabilities and the integers, you get the expected value of the distribution. So, if you have a mean value in mind, say $k$, you can solve the following equation:
$$k = 1\times P(1) + 2\times P(2) + 3\times P(3) + 4\times P(4)$$
You can arbitrarily choose two of the probabilities and solve for the third, which determines the fourth (because $P(1)=1-(P(2)+P(3)+P(4))$ because the probabilities must sum to $1$). For example, let $k=2.3$, $P(4)=.1$, and $P(3)=.2$. Then we have that 
$$k = 1 \times [1-(P(2)+P(3)+P(4)] + 2\times P(2) + 3\times P(3) + 4\times P(4)$$
$$2.3 = [1 - (P(2)+.1+.2)] + 2*P(2) + 3\times .2 + 4\times .1$$
$$2.3 = .7 + P(2) + .6 + .4$$
$$P(2)=.6$$
$$P(1)=1-(P(2)+P(3)+P(4)=1 - (.6+.1+.2)=.1$$
So you can run x <- sample(c(1, 2, 3, 4), 1e6, replace = TRUE, prob = c(.1, .6, .2, .1)) and mean(x) is approximately $2.3$
A: The question is under-specified in that the constraints on the frequencies
\begin{align}n_1+2n_2+3n_3+4n_4&=100M\\n_1+n_2+n_3+n_4&=100\end{align}
do not determine a distribution: "random" is not associated with a particular distribution, unless the OP means "uniform". For instance, if there exists one solution $(n_1^0,n_2^0,n_3^0,n_4^0)$ to the above system, then the distribution degenerated at this solution is producing a random draw that is always $(n_1^0,n_2^0,n_3^0,n_4^0)$.
In the case the question is about simulating a Uniform distribution over the grid\begin{align}n_1+2n_2+3n_3+4n_4&=100M\\n_1+n_2+n_3+n_4&=100\end{align}one can always use a Metropolis-Hastings algorithm. Starting from $(n_1^0,n_2^0,n_3^0,n_4^0)$, create a Markov chain by proposing symmetric random perturbations of the vector $(n_1^t,n_2^t,n_3^t,n_4^t)$ and accept if the result is within $\{1,2,3,4\}^4$ and satisfies the constraints.
For instance, here is a crude R rendering:
cenM=293
#starting point (n¹,n³,n⁴)
n<-sample(1:100,3,rep=TRUE)
while((sum(n)>100)|(n[2]-n[1]+2*n[3]!=cenM-200))
    n<-sample(1:100,3,rep=TRUE)
#Markov chain
for (t in 1:1e6){
  prop<-n+sample(-10:10,3,rep=TRUE)
  if ((sum(prop)<101)&
      (prop[2]-prop[1]+2*prop[3]==cenM-200)&
      (min(prop)>0)) 
        n=prop}
c(n[1],100-sum(n),n[-1])

with the distribution of $(n_1,n_3,n_4)$ over the 10⁶ iterations:

In case you want draws of the integers themselves,
 sample(c(rep(1,n[1]),rep(2,100-sum(n)),rep(3,n[2]),rep(4,n[3])))

is a quick & dirty way to produce a sample.
A: Here is a simple algorithm: Create $n-1$ random integers in the range $[1,4]$ and calculate the $n^{th}$ integer for the mean to be equal to the specified value. If that number is smaller than $1$ or larger than $4$, then one by one distribute the surplus/lacking onto other integers, e.g. if the integer is $5$, we have $1$ surplus; and we may add this to the next integer if it's not $4$, else add to the next etc. Then, shuffle the entire array.
A: As a supplement to whuber's answer, I've written a script in Python which goes through each step of the sampling scheme. Note that this is meant for illustrative purposes and is not necessarily performant.
Example output:
n=10, s=20, k=4

Starting grid
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
X X X X X X X X X X

Filled in grid
X X . . X . X . . X
. . X X X . . . . .
. . . . X X . . . .
X X X X X X X X X X

Final grid
X X . . X . X . . X
. . X X X . . . . .
. . . . X X . . . .
X X X X X X X X X X
2 2 2 2 4 2 2 1 1 2

The script:
import numpy as np

# Define the starting parameters
integers = [1, 2, 3, 4]
n = 10
s = 20
k = len(integers)


def print_grid(grid, title):
    print(f'\n{title}')
    for row in grid:
        print(' '.join([str(element) for element in row]))


# Create the starting grid
grid = []
for i in range(1, k + 1):
    if i < k:
        grid.append(['.' for j in range(n)])
    else:
        grid.append(['X' for j in range(n)])

# Print the starting grid
print_grid(grid, 'Starting grid')

# Randomly and uniformly fill in the remaining rows
indexes = np.random.choice(range((k - 1) * n), s - n, replace=False)
for i in indexes:
    row = i // n
    col = i % n
    grid[row][col] = 'X'

# Print the filled in grid
print_grid(grid, 'Filled in grid')

# Compute how many cells were selected in each column
column_counts = []
for col in range(n):
    count = sum(1 for i in range(k) if grid[i][col] == 'X')
    column_counts.append(count)
grid.append(column_counts)

# Print the final grid and check that the column counts sum to s
print_grid(grid, 'Final grid')
print()
print(f'Do the column counts sum to {s}? {sum(column_counts) == s}.')

A: I've turned whuber's answer into an r function.  I hope it helps someone.


*

*n is how many integers you want;

*t is the mean you want; and 

*k is the upper limit you want for your returned values


whubernator<-function(n=NULL, t=NULL, kMax=5){
  z = tabulate(sample.int(kMax*(n), (n)*(t),replace =F) %% (n)+1, (n))
  return(z)
}

It seems to work as expected:
> w = whubernator(n=10,t=4.2)
> mean(w)
[1] 4.2
> length(w)
[1] 10
> w
 [1] 3 5 3 5 5 3 4 5 5 4

It can return 0s, which matches my needs.
> whubernator(n=2,t=0.5)
[1] 1 0

