Sono d'accordo con X'ian sul fatto che il problema sia sotto specificato. Tuttavia, esiste una soluzione elegante, scalabile, efficiente, efficace e versatile che vale la pena considerare.

Poiché il prodotto della media e della dimensione del campione è uguale alla somma del campione, il problema riguarda la generazione di un campione casuale di valori $ n $ nell'insieme $ \ {1,2, \ ldots, k \} $ che somma a $ s $ (assumendo $ n \ le s \ le kn, $ ovviamente).

Per spiegare la soluzione proposta e, spero, giustificare l'affermazione di elegance, offro un'interpretazione grafica di questo schema di campionamento. Disponi una griglia di $ k $ righe e $ n $ colonne. Seleziona ogni cella nella prima riga. Seleziona in modo casuale (e uniforme) $ sn $ delle celle rimanenti nelle righe da $ 2 $ a $ k. $ Il valore dell'osservazione $ i $ nel campione è il numero di celle selezionate nella colonna $ i: $

Questa griglia $ 4 \ times 100 $ è rappresentata da punti neri nelle celle non selezionate e da toppe colorate nelle celle selezionate. È stato generato per produrre un valore medio di $ 2, $ quindi $ s = 200. $ Quindi, $ 200-100 = 100 $ sono state selezionate casualmente tra le prime $ k-1 = 3 $ righe. I colori rappresentano i numeri delle celle selezionate in ogni colonna. Ci sono $ 28 $ , $ 47 $ due, $ 22 $ tre e $ 3 $ quattro. Il campione ordinato corrisponde alla sequenza di colori dalla colonna $ 1 $ alla colonna $ n = 100. $

Per dimostrare la scalabilità e l'efficienza, ecco un comando R per generare un campione secondo questo schema. La domanda riguarda il caso in cui $ k = 4, n = 100 $ e $ s $ è $ n $ volte la media desiderata del campione:

  tabulate (sample.int ((k-1) * n, s-n) %% n + 1, n) + 1
 

Perché sample.int richiede $ O (sn) $ tempo e $ O ( (k-1) n) $ spazio e tabulate richiede $ O (n) $ tempo e spazio, questo algoritmo richiede $ O (\ max (sn, n)) $ tempo e $ O (kn) $ spazio : questo è scalable. Con $ k = 4 $ e $ n = 100 $ la mia workstation impiega solo 12 microsecondi per eseguire questo calcolo : questo è efficient.

(Ecco una breve spiegazione del codice. Tieni presente che i numeri interi $ x $ in $ \ {1,2, \ ldots, (k-1) n \} $ può essere espresso in modo univoco come $ x = nj + i $ dove $ j \ in \ {0,1, \ ldots, k-2 \} $ e $ i \ in \ {1,2, \ ldots, n \}. $ Il codice prende un esempio di tali $ x, $ li converte nel loro $ ( i, j) $ coordinate della griglia, conta quante volte ogni $ i $ appare (che sarà compreso tra $ 0 $ a $ k-1 $ ) e aggiunge $ 1 $ a ogni conteggio.)

Perché questo può essere considerato efficace? Una ragione è che le proprietà distributive di questo schema di campionamento sono semplici da elaborare:

• È scambiabile: tutte le permutazioni di qualsiasi campione sono ugualmente probabili.

• La possibilità che il valore $ x \ in \ {1,2, \ ldots, k \} $ appaia nella posizione $ i, $ che scriverò come $ \ pi_i (x), $ è ottenuto tramite un argomento di conteggio ipergeometrico di base come $$ \ pi_i (x) = \ frac {\ binom {k-1} {x-1} \ binom {(n-1) (k-1)} {sn-x +1}} {\ binom {n (k-1)} {sn}}. $$ Ad esempio, con $ k = 4, $ $ n = 100, $ e una media di $ 2,0 $ (in modo che $ s = 200 $ ) è probabile che $ \ pi = (0.2948, 0.4467, 0.2222, 0.03630), $ concordino strettamente con le frequenze nel campione precedente. Ecco i grafici di $ \ pi_1 (1), \ pi_1 (2), \ pi_1 (3), $ e $ \ pi_1 (4) $ in funzione della somma:

  

• La possibilità che il valore $ x $ appaia nella posizione $ i $ mentre il il valore $ y $ appare nella posizione $ j $ si trova in modo simile a $$ \ pi_ {ij} (x, y) = \ frac {\ binom {k-1} {x-1} \ binom {k-1} {y-1} \ binom {(n-1) (k-1)} {snx-y + 2}} {\ binom {n (k-1)} {sn}}. $$

Queste probabilità $ \ pi_i $ e $ \ pi_ {ij} $ consentono di applicare il stimatore Horvitz-Thompson per questo progetto di campionamento probabilistico e per calcolare i primi due momenti delle distribuzioni di varie statistiche.

Infine, questa soluzione è versatile nella misura in cui consente variazioni semplici e prontamente analizzabili per controllare la distribuzione del campionamento.Ad esempio, è possibile selezionare celle sulla griglia con probabilità specificate ma disuguali in ogni riga o con un modello simile a un'urna per modificare le probabilità man mano che il campionamento procede, controllando così le frequenze dei conteggi delle colonne.

La domanda è sotto specificata in quanto i vincoli sulle frequenze \ begin {align} n_1 + 2n_2 + 3n_3 + 4n_4& = 100M \\ n_1 + n_2 + n_3 + n_4& = 100 \ end {align} non determinare una distribuzione: "casuale" non è associato a una particolare distribuzione, a meno che l'OP non significhi "uniforme". Ad esempio, se esiste una soluzione $ (n_1 ^ 0, n_2 ^ 0, n_3 ^ 0, n_4 ^ 0) $ per il sistema precedente, la distribuzione è degenerata a questa soluzione si produce un'estrazione casuale che è sempre $ (n_1 ^ 0, n_2 ^ 0, n_3 ^ 0, n_4 ^ 0) $ .

Nel caso in cui la domanda riguardi la simulazione di una distribuzione uniforme sulla griglia \ begin {align} n_1 + 2n_2 + 3n_3 + 4n_4& = 100M \\ n_1 + n_2 + n_3 + n_4& = 100 \ end {align} si può sempre usare un algoritmo Metropolis-Hastings. A partire da $ (n_1 ^ 0, n_2 ^ 0, n_3 ^ 0, n_4 ^ 0) $ , crea una catena di Markov proponendo perturbazioni casuali simmetriche del vettore $ (n_1 ^ t, n_2 ^ t, n_3 ^ t, n_4 ^ t) $ e accetta se il risultato è compreso tra $ \ {1,2,3,4 \} ^ 4 $ e soddisfa i vincoli.

Ad esempio, ecco un semplice rendering R:

  cenM = 293
# punto di partenza (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 catena
for (t in 1: 1e6) {
  prop<-n + campione (-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-somma (n), n [-1])
 

con la distribuzione di $ (n_1, n_3, n_4) $ nelle iterazioni 10⁶:

Nel caso in cui desideri estrarre gli interi stessi,

  sample (c (rep (1, n [1]), rep (2.100-sum (n)), rep (3, n [2]), rep (4, n [3])) )
 

è un modo rapido e sporco di & per produrre un campione.

Voglio ... uh ... "attenuare" la straordinaria risposta di @ whuber, che @TomZinger dice è troppo difficile da seguire. Con questo voglio dire che voglio ridescriverlo in termini che penso che Tom Zinger capirà, perché è chiaramente la migliore risposta qui. E mentre Tom usa gradualmente il metodo e scopre che ha bisogno, diciamo, di conoscere la distribuzione dei campioni piuttosto che solo la loro media, la risposta di chiunque sarà proprio quello che sta cercando.

In breve: non ci sono idee originali qui, solo una spiegazione più semplice.

Desideri creare $ n $ numeri interi da $ 1 $ a $ 4 $ con media $ r $ . Suggerisco di calcolare $ n $ numeri interi da $ 0 $ a $ 3 $ con mean $ r-1 $ , quindi aggiungendone uno a ciascuno di essi. Se puoi fare quest'ultima cosa, puoi risolvere il primo problema. Ad esempio, se vogliamo 10 numeri interi compresi tra $ 1 $ e $ 4 $ con media $ 2,6 $ , possiamo annotare questi $ 10 $ interi tra $ 0 $ e $ 3 $ ...

0,3,2,1,3,1,2,1,3,0

la cui media è $ 1,6 $ ; se aumentiamo ciascuno di $ 1 $ , otteniamo

1,4,3,2,4,2,3,2,4,1

la cui media è $ 2,6 $ . È così semplice.

Ora pensiamo ai numeri da $ 0 $ a $ 3 $ . Li considererò come "quanti oggetti ho in un set" piccolo "?" Potrei non avere elementi, un elemento, due elementi o tre elementi. Quindi l'elenco

0,3,2,1,3,1,2,1,3,0

rappresenta dieci diversi piccoli set. Il primo è vuoto; il secondo ha tre elementi e così via. Il numero totale di elementi in tutti gli insiemi è la somma dei dieci numeri, ovvero $ 16 $ . E il numero medio di elementi in ogni set è questo totale, diviso per $ 10 $ , quindi $ 1,6 $ .

L'idea di chiunque sia questa: supponi di farti dieci piccoli set, con il numero totale di elementi $ 10t $ per un certo numero $ t $ . Quindi la dimensione media degli insiemi sarà esattamente $ t $ . Allo stesso modo, se imposti $ n $ set con un numero totale di elementi $ nt $ span >, il numero medio di elementi in un set sarà $ t $ . Dici di essere interessato al caso $ n = 100 $ .

Rendiamolo concreto per il tuo esempio: vuoi 100 elementi tra 1 e 4 la cui media è $ 1,9 $ . Usando l'idea del mio primo paragrafo, lo cambierò in "make $ 100 $ ints tra $ 0 $ span> e $ 3 $ la cui media è $ 0,9 $ ". Quando ho finito, aggiungerò $ 1 $ a ciascuno dei miei int per ottenere una soluzione al tuo problema. Quindi la mia media target è $ t = 0.9 $ .

Voglio creare set di $ 100 $ , ciascuno con un valore compreso tra $ 0 $ e $ 3 $ elementi in esso, con una dimensione media del set di $ 0,9 $ .

Come ho osservato sopra, ciò significa che deve esserci un totale di $ 100 \ cdot 0.9 = 90 $ elementi in i set. Dai numeri $ 1, 2, \ ldots, 300 $ , selezionerò esattamente $ 90 $ . Posso indicare quelli selezionati facendo un elenco di 300 punti e X:

..X .... X ... XX ...

dove l'elenco sopra indica che ho selezionato i numeri 3, 9, 13, 14 e poi molti altri che non ho mostrato perché mi sono stufato di digitare. :) Posso prendere questa sequenza di 300 punti e X e suddividerla in tre gruppi di 100 punti ciascuno, che dispongo uno sopra l'altro, ottenendo qualcosa di simile a questo:

  ... X .... X..X ..... X ...
.X ... X ..... X ... X .....
..X ... X.X..X ...... X ..
 

ma continua per 100 articoli completi in ogni riga. Il numero di X in ogni riga potrebbe essere diverso: potrebbero esserci 35 nella prima riga, 24 nella seconda e 31 nella terza, ad esempio, e va bene. [Grazie a Whuber per aver sottolineato che avevo sbagliato in una prima bozza!]

Ora guarda ogni colonna : ogni colonna può essere considerata come un insieme e quell'insieme contiene da 0 a 3 "X". Posso scrivere i conteggi sotto le righe per ottenere qualcosa del genere:

  ... X .... X..X ..... X ...
.X ... X ..... X ... X .....
..X ... X.X..X ...... X ..
011101102003000101100
 

Vale a dire, ho prodotto 100 numeri, ciascuno compreso tra 1 e 3. E la somma di questi 100 numeri deve essere il numero di X, totale, in tutte e tre le righe, che era 90. Quindi la media deve essere $ 90/100 = 0,9 $ , come desiderato.

Ecco quindi i passaggi per ottenere 100 numeri interi compresi tra 1 e 4 la cui media è esattamente $ s $ .

1. Lascia $ t = s - 1 $ .
2. Calcola $ k = 100 t $ ; è il numero di X che inseriremo nelle righe, in totale.
3. Crea un elenco di 300 punti o X, di cui $ k $ .
4. Dividilo in tre righe di 100 punti o X, ciascuna contenente circa un terzo delle X, più o meno.
5. Disporli in un array e calcolare le somme delle colonne, ottenendo 100 numeri interi compresi tra $ 0 $ e $ 3 $ span >. La loro media sarà $ t $ .
6. Aggiungine uno a ogni somma di colonna per ottenere 100 numeri interi compresi tra $ 1 $ e $ 4 $ la cui media è $ s $ .

Ora la parte complicata di questo è proprio nel passaggio 4: come si selezionano $ 300 $ articoli, $ k $ di cui sono "X" e gli altri $ 300-k $ di cui sono "."? Bene, risulta che R ha una funzione che fa esattamente questo.

E poi whuber ti dice come usarlo: scrivi tu

  tabulate (sample.int ((k-1) * n, s-n) %% n + 1, n)
 

Nel tuo caso particolare, $ n = 100 $ e $ s $ , il numero totale di elementi in tutti i piccoli insiemi, è $ 100r $ e desideri numeri compresi tra $ 1 $ e $ 4 $ , quindi $ k = 4 $ , quindi $ k -1 $ (la dimensione massima per un 'insieme piccolo') è 3, quindi diventa

  tabulate (sample.int (3 * 100, 100r-100) %% 100 + 1, n)
 

o

  tabulate (sample.int (3 * 100, 100 * (r-1)) %% 100 + 1, 100)
 

oppure, utilizzando il mio nome $ t $ per $ r - 1 $ , diventa

  tabulate (sample.int (3 * 100, 100 * t) %% 100 + 1, 100)
 

Il "+1" alla fine della sua formula originale è esattamente il passaggio necessario per convertire da "numeri compresi tra $ 0 $ e $ 3 $ " a "numeri compresi tra $ 1 $ e $ 4 $ ".

Lavoriamo dall'interno verso l'esterno e semplifichiamo in $ n = 10 $ in modo da poter mostrare output di esempio:

  tabulate (sample.int (3 * 10, 10 * t) %% 10 + 1, 10)
 

E miriamo a $ t = 1.9 $ , quindi questo diventa

  tabulate (sample.int (3 * 10, 10 * 1.9) %% 10 + 1, 10)
 

A partire da sample.int (3 * 10, 10 * 1.9) : questo produce un elenco di $ 19 $ interi compresi tra $ 1 $ e $ 30 $ . (cioè, ha risolto il problema di scegliere $ k $ numeri dal tuo totale - $ 300 $ in il tuo vero problema, $ 30 $ nel mio esempio più piccolo).

Come ricorderai, vogliamo produrre tre righe di dieci punti e X ciascuna, qualcosa come

  X.X.XX.XX.
 XXXX.XXX ..
 XX.X.XXX ..
 

Possiamo leggere questo da sinistra a destra dall'alto verso il basso (cioè, normale ordine di lettura) per produrre un elenco di posizioni per X: il primo elemento è un punto; la seconda e la terza sono X e così via, quindi il nostro elenco di posizioni inizia con $ 1, 3, 5, 6, \ ldots $ . Quando arriviamo alla fine di una riga, continuiamo a contare, quindi per l'immagine sopra, le posizioni X sarebbero $ 1, 3, 5, 6, 8, 9, 11, 12, 13, 14, 16, 17, 18, 21, 22, 24, 26, 27, 28 $ . È chiaro?

Ebbene, il codice Whubers produce esattamente quell'elenco di posizioni con la sua sezione più interna.

L'elemento successivo è %% 10 ; che prende un numero e produce il suo resto sulla divisione per dieci. Quindi il nostro elenco diventa $ 1, 3, 5, 6, 8, 9, 1, 2, 3, 4, 6, 7, 8, 1, 2, 4, 6, 7, 8 $ . Se lo suddividiamo in tre gruppi --- quelli che provengono da numeri compresi tra $ 1 $ e $ 10 $ , quelli che provenivano da numeri da $ 11 $ a $ 20 $ e quelli che provenivano da numeri $ 21 $ a $ 30 $ , otteniamo $ 1, 3, 5, 6, 8, 9 $ , quindi $ 1, 2, 3, 4, 6, 7, 8, $ e infine $ 1, 2, 4, 6, 7, 8 $ . Quelli ti dicono dove sono le X in ciascuna delle tre righe. C'è un sottile problema qui: se ci fosse stata una X nella posizione 10 nella prima riga, il primo dei nostri tre elenchi sarebbe stato $ 1, 3, 5, 6, 8, 9 , 0 $ e alla funzione tabulate non piace "0". Quindi whuber aggiunge 1 a ogni elemento nell'elenco per ottenere $ 2, 4, 6, 7, 9, 10, 1 $ . Passiamo al calcolo complessivo:

  tabulate (sample.int (3 * 10, 10 * 1.9) %% 10 + 1, 10)
 

Questo richiede "quei numeri $ 30 $ , ciascuno che indica se c'è una X in una colonna, dimmi quante volte ogni colonna (da $ 1 $ a $ 10 $ --- questo è quello che ti dice il" 10 "finale) appare, cioè dimmi quante X sono in ogni colonna. Il risultato è 0 3 2 2 2 1 3 2 3 1 che (a causa dello spostamento di una cosa) devi leggere come "non ci sono X nella decima colonna; ci sono 3 X nella prima colonna; ci sono 2 X nella seconda colonna" e così via su a "c'è una X nella nona colonna".

Ciò fornisce dieci numeri interi compresi tra $ 0 $ e $ 3 $ la cui somma è $ 19 $ , quindi la cui media è $ 1,9 $ . Se aumenti ciascuno di 1, ottieni dieci numeri interi compresi tra $ 1 $ e $ 4 $ la cui somma è $ 29 $ , quindi un valore medio di $ 2,9 $ .

Puoi generalizzare a $ n = 100 $ , spero.

Puoi utilizzare sample () e selezionare probabilità specifiche per ogni numero intero. Se si somma il prodotto delle probabilità e degli interi, si ottiene il valore atteso della distribuzione. Quindi, se hai in mente un valore medio, ad esempio $ k $ , puoi risolvere la seguente equazione: $$ k = 1 \ times P (1) + 2 \ times P (2) + 3 \ times P (3) + 4 \ times P (4) $$ span > Puoi scegliere arbitrariamente due delle probabilità e risolvere la terza, che determina la quarta (perché $ P (1) = 1- (P (2) + P (3) + P (4)) $ perché le probabilità devono sommarsi a $ 1 $ ). Ad esempio, lascia $ k = 2.3 $ , $ P (4) =. 1 $ e $ P (3) =. 2 $ . Allora abbiamo quello $$ k = 1 \ times [1- (P (2) + P (3) + P (4)] + 2 \ times P (2) + 3 \ times P ( 3) + 4 \ volte 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 $$

Quindi puoi eseguire x <- sample (c (1, 2, 3, 4), 1e6, replace = TRUE, prob = c (.1, .6, .2, .1)) e mean (x) è approssimativamente $ 2,3 $

Ecco un semplice algoritmo: crea $ n-1 $ numeri interi casuali nell'intervallo $ [1,4]$ e calcola l'intero $ n ^ {th} $ affinché la media sia uguale al valore specificato.Se il numero è inferiore a $ 1 $ o superiore a $ 4 $ , distribuire uno per uno il surplus /privo di altri numeri interi, ad esse il numero intero è $ 5 $ , abbiamo $ 1 $ surplus;e possiamo aggiungerlo al numero intero successivo se non è $ 4 $ , altrimenti aggiungerlo al successivo ecc. Quindi, mescola l'intero array.

Come supplemento alla risposta di chiunque, ho scritto uno script in Python che segue ogni passaggio dello schema di campionamento. Nota che questo è inteso a scopo illustrativo e non è necessariamente performante.

Risultato di esempio:

  n = 10, s = 20, k = 4

Griglia di partenza
. . . . . . . . . .
. . . . . . . . . .
. . . . . . . . . .
X X X X X X X X X X

Griglia riempita
X X. . X . X . . X
. . X X X. . . . .
. . . . X X. . . .
X X X X X X X X X X

Griglia finale
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
 

Lo script:

  importa numpy come np

# Definisce i parametri di partenza
numeri interi = [1, 2, 3, 4]
n = 10
s = 20
k = len (numeri interi)


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


# Crea la griglia di partenza
griglia = []
per i nell'intervallo (1, k + 1):
    se io < k:
        grid.append (['.' for j in range (n)])
    altro:
        grid.append (['X' per j nell'intervallo (n)])

# Stampa la griglia di partenza
print_grid (griglia, 'Griglia iniziale')

# Compila in modo casuale e uniforme le righe rimanenti
indexes = np.random.choice (range ((k - 1) * n), s - n, replace = False)
per i negli indici:
    riga = i // n
    col = i% n
    griglia [riga] [col] = "X"

# Stampa la griglia compilata
print_grid (grid, 'Filled in grid')

# Calcola quante celle sono state selezionate in ogni colonna
column_counts = []
per col nell'intervallo (n):
    count = sum (1 for i in range (k) if grid [i] [col] == 'X')
    column_counts.append (count)
grid.append (column_counts)

# Stampa la griglia finale e controlla che la colonna conti la somma per s
print_grid (griglia, 'Griglia finale')
Stampa()
print (f'Do la colonna conta somma a {s}? {sum (column_counts) == s}. ')
 

Ho trasformato la risposta di whuber in una funzione r.Spero che aiuti qualcuno.

• n è il numero di numeri interi che desideri;
• t è la media che vuoi;e
• k è il limite massimo che desideri per i valori restituiti

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

Sembra funzionare come previsto:

  > w = whubernator (n = 10, t = 4.2)
> mean (w)
[1] 4.2
> lunghezza (w)
[1] 10
> w
 [1] 3 5 3 5 5 3 4 5 5 4
 

Può restituire 0, che corrisponde alle mie esigenze.

  > whubernator (n = 2, t = 0,5)
[1] 1 0