Simulation of IQ imbalance using MATLAB
Bottom line:

You cannot uniquely determine $\alpha$ and $\theta$ from IRR given in dB since it has a dependence on both given by the following approximation:

$$IRR \approx 20Log_{10}\bigg(\bigg|\frac{1-  \alpha - j\theta}{1+  \alpha + j\theta}\bigg|\bigg)$$

Where

$IRR$ is the image rejection ratio in dB

$\alpha$ is amplitude imbalance as a magnitude ratio (typically 0.9 to 1.1)

$\theta$ is phase imbalance in radians (typically -0.1 to 0.1)

This uses the small angle approximation such that $sin(\theta) \approx \theta$ which covers most IRR applications. 

For cases where there is no quadature imbalance and only amplitude imbalance the following relationship is exact and you can then solve for $\alpha$ from IRR:

$$IRR = 20Log_{10}\bigg(\bigg|\frac{1-  \alpha}{1+  \alpha}\bigg|\bigg)$$

From which:

$$\alpha = 1\pm \bigg(1-\frac{1-r}{1+r}\bigg)$$

Where 
$r = 10^{IRR/20}$

With IRR given as a negative quantity indicating the image rejection ratio in dB.

For example the amplitude imbalance to provide -30  dB image rejection would be as follows assuming perfect quadrature:

$r = 10^{-30/20} = 0.0316$

$\alpha = 0.9387, 1.0613$

(Or +.51 dB, -.54 dB)

A similar approximation can be given for small angle $\theta$ for when there is no amplitude imbalance. Bessel functions can be used to derive an exact result for $\theta$ applicable to larger angles.

Details

If you complex waveform is already generated (as it would be trivial to generate otherwise) you can use the corrections as given in the diagram below (if it corrects a phase and amplitude error it will equivalently induce a phase and amplitude error):

So multiply I by your $\theta$ term and add to Q to add quadrature error ($\theta$ will be approximately your angle in radians using small angle approximations), and scale I by your $\alpha$ term to add amplitude error. 

$$I_{out} = \alpha I_{in}$$
$$Q_{out} = \theta I_{in} + Q_{in}$$



Both $\alpha$ and $\theta$ will both effect IRR. If the IRR is given in dB, there is no way to uniquely determine $\alpha$ and $\theta$. Note the following to see how IRR can be determined from knowing $\alpha$ and $\theta$:

Graphically it is easy to derive this since a circle represents perfect I and Q balance and results in image rejection. We can see this simply by using Euler's Identities:

Given Euler's Identities as

$$cos(\omega t) = \frac{e^{j\omega t}}{2} + \frac{e^{-j\omega t}}{2}$$
$$sin(\omega t) = \frac{e^{j\omega t}}{2j} - \frac{e^{-j\omega t}}{2j}$$

And a circle with perfect IQ balance as:
$$y(\theta) = cos(\omega t) + jsin(\omega t)$$

$$= \frac{e^{j\omega t}}{2} + \frac{e^{-j\omega t}}{2} +  j\frac{e^{j\omega t}}{2j} - j\frac{e^{-j\omega t}}{2j}$$

$$= \frac{e^{j\omega t}}{2} + \frac{e^{-j\omega t}}{2} +  \frac{e^{j\omega t}}{2} - \frac{e^{-j\omega t}}{2} = e^{j\omega t}$$

Resulting in the well known Euler's relationship as 

$$ e^{j\omega t} = cos(\omega t) + jsin(\omega t)$$

The above is insightful as it demonstrates the sideband rejection through phase summation from the I and Q components being in perfect balance and quadrature, and recognizing that $e^{j\omega t}$ represents a single sideband tone (a positive frequency $\omega$ while a sinusoid has two sidebands or a positive and negative frequency). 

So similarly we can evaluate the effect of amplitude and phase imbalance. From the block diagram above showing the relationship of $\alpha$ and $\theta$ for small angle approximation when $sin(\theta) \approx \theta$:

$$y(\theta) = \alpha cos(\omega t) + jsin(\omega t) + j\theta cos(\omega t)$$

$$y(\theta) = \frac{ \alpha e^{j\omega t}}{2} + \frac{  \alpha e^{-j\omega t}}{2} + \frac{je^{j\omega t}}{2j} - \frac{je^{-j\omega t}}{2j} + \frac{j \theta e^{j\omega t}}{2} + \frac{ j \theta e^{-j\omega t}}{2} $$

$$y(\theta) = \frac{ \alpha e^{j\omega t}}{2} + \frac{  \alpha e^{-j\omega t}}{2} + \frac{e^{j\omega t}}{2} - \frac{e^{-j\omega t}}{2} + \frac{j \theta e^{j\omega t}}{2} + \frac{ j \theta e^{-j\omega t}}{2} $$

$$y(\theta) = \frac{(1+ \alpha +j\theta)}{2} e^{j\omega t} - \frac{(1-  \alpha - j\theta)}{2} e^{-j\omega t} $$

The IRR in dB would be $20Log_{10}$ of the ratio of the two exponential coefficients given as:

$$IRR \approx 20Log_{10}\bigg(\bigg|\frac{1-  \alpha - j\theta}{1+  \alpha + j\theta}\bigg|\bigg)$$