<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def theo1(data, rate=1.0, data_type="phase", taus=None):
 """ PRELIMINARY - REQUIRES FURTHER TESTING. Theo1 is a two-sample variance with improved confidence and extended averaging factor range. .. math:: \\sigma^2_{THEO1}(m\\tau_0) = { 1 \\over (m \\tau_0 )^2 (N-m) } \\sum_{i=1}^{N-m} \\sum_{\\delta=0}^{m/2-1} {1\\over m/2-\\delta}\\lbrace ({x}_{i} - x_{i-\\delta +m/2}) + (x_{i+m}- x_{i+\\delta +m/2}) \\rbrace^2 Where :math:`10<=m<=N-1` is even. FIXME: bias correction NIST SP 1065 eq (30) page 29 Parameters data: np.array Input data. Provide either phase or frequency (fractional, adimensional). rate: float The sampling rate for data, in Hz. Defaults to 1.0 data_type: {'phase', 'freq'} Data type, i.e. phase or frequency. Defaults to "phase". taus: np.array Array of tau values, in seconds, for which to compute statistic. Optionally set taus=["all"|"octave"|"decade"] for automatic tau-list generation. """

    phase = input_to_phase(data, rate, data_type)

    tau0 = 1.0/rate
    (phase, ms, taus_used) = tau_generator(phase, rate, taus, even=True)

    devs = np.zeros_like(taus_used)
    deverrs = np.zeros_like(taus_used)
    ns = np.zeros_like(taus_used)

    N = len(phase)
    for idx, m in enumerate(ms):
        m = int(m) # to avoid: VisibleDeprecationWarning: using a
                   # non-integer number instead of an integer will
                   # result in an error in the future
        assert m % 2 == 0 # m must be even
        dev = 0
        n = 0
        for i in range(int(N-m)):
            s = 0
            for d in range(int(m/2)): # inner sum
                pre = 1.0 / (float(m)/2 - float(d))
                s += pre*pow(phase[i]-phase[i-d+int(m/2)] +
                             phase[i+m]-phase[i+d+int(m/2)], 2)
                n = n+1
            dev += s
        assert n == (N-m)*m/2 # N-m outer sums, m/2 inner sums
        dev = dev/(0.75*(N-m)*pow(m*tau0, 2))
        # factor 0.75 used here? http://tf.nist.gov/general/pdf/1990.pdf
        # but not here? http://tf.nist.gov/timefreq/general/pdf/2220.pdf page 29
        devs[idx] = np.sqrt(dev)
        deverrs[idx] = devs[idx] / np.sqrt(N-m)
        ns[idx] = n

    return remove_small_ns(taus_used, devs, deverrs, ns)