David P Tew

Wim Klopper
email: willem.klopper@kit.edu

Open-shell explicitly correlated F12 methods

Keywords: Open-shell systems, Electron correlation, Geminals, Correlation cusp

In calculations on open-shell systems, spin-flipped geminals must be used in explicitly correlated F12 methods with spin-unrestricted or restricted open-shell Hartree-Fock reference wave functions.

We suggest to use the sp ansatz of Bokhan et al. [J. Chem. Phys. 131, 084105 (2009)] for calculations with fixed F12 amplitudes, but to use a novel contracted geminal approach if the F12 amplitudes are to be optimised. This new approach is denoted MP2-F12-o respectively CCSD(F12)-o. The performance of the new approach is assessed by calculating the atomisation energies of a test set comprising 106 molecules containing the atoms H, C, N, O and F, and by calculating the ionisation energies and electron affinities of the atoms C, N, O and F.

R12 explicitly-correlated electronic structure methods overcome the slow convergence of orbital expansions by employing geminal basis functions that closely resemble the correlation holes in electronic wave functions. 1,2 These methods have developed rapidly over the last seven years and are emerging as robust tools for computational chemistry, where small basis sets may be used to obtain near basis set limit accuracy. 3 The principal differences between the early and modern R12 methods are the use of a complementary auxiliary basis set (CABS) for the resolution of the identity approximation, employed to evaluate the manyelectron integrals, 4,5 and the use of a Slater-type correlation factor for the geminal basis functions. 6 Methods that use non-linear correlation factors instead of the original linear factor are commonly referred to as F12 methods. So far, the F12 approach has been combined with MP2, 6 CASPT2, 7 and coupled-cluster wave function methods up to CCSDTQ. [8][9][10][11] In addition, many researchers have improved the efficiency of F12 methods through density fitting, 12,13 improved methods for evaluating Fock matrix elements for geminal functions, 14,15 localisation techniques, [16][17][18][19] and by constructing specially optimised basis sets. [20][21][22] One further development is the rational generator approach of Ten-no. 23 In this approach the s-and p-wave coalescence conditions 24,25 for the first-order wave function are used to predetermine the F12 amplitudes instead of optimising them. The fixed-amplitude method is computationally less demanding because only the diagonal elements of the spin-adapted F12 matrix elements need to be computed and no equations for the F12 amplitudes need to be solved. Moreover, the fixed-amplitude method is free from geminal basis set superposition error, which arises if the amplitudes are optimised, 26 and is free from the numerical problems associated with inverting the geminal Fock matrix, which becomes non-positive definite in some limits. 27 The loss of accuracy due to the reduced variational degrees of freedom is slight and the fixed-amplitude method is used increasingly in F12 calculations, particularly for explicitly-correlated coupled-cluster singles and doubles methods (CCSD-F12). [28][29][30][31][32][33][34] The extension of the fixed-amplitude approach to open-shell systems has been attempted only recently 27,30,32,35 and has led to a re-evaluation of the geminal basis used for open-shell calculations. 32 Bokhan et al. demonstrated that, contrary to the closed-shell case, it is not possible to satisfy both the s-and p-wave UMP2 coalescence conditions using the functions f (r 12 )|ij alone. |ij is the usual two-electron determinant of occupied spin-orbitals i and j have subsequently performed calculations at the UMP2-F12 level. 35 They demonstrate that the basis set convergence without the spin-flipped basis is limited by the X -5 behaviour for the pair energies where one orbital is singly occupied in the HF reference, arising from the unsatisfied p-wave derivative discontinuity. By including the spin-flipped basis, the X -7 convergence behaviour of F12 methods for closed shell species is achieved for openshell species. For high spin ROHF references, Knizia et al. have also included spin-flipped geminal basis functions in their RMP2-F12 and CCSD-F12x methods, 30 noting improved performance for atomisation energies, electron affinities and ionisation potentials compared

to their earlier open-shell fixed-amplitude approach. 27 Nevertheless, a number of questions remain unanswered for open-shell F12 calculations.

For UHF references, the combined basis of normal and spin-flipped geminals contains many near linear-dependencies and including both sets with full variational flexibility leads to severe numerical problems. Indeed, Bokhan et al. only considered the diagonal orbital-variant method and, even here, a singular value decomposition was necessary. We consider that optimising coefficients for both the normal and spin-flipped geminals is undesirable since it results in an unbalanced treatment of open-and closed-shell systems, favouring the openshell species. It is also numerically problematic and can lead to steps in potential energy surfaces. In this work we suggest that a contracted geminal basis should be used, that is, fixed linear combinations of normal and spin-flipped geminals with contraction coefficients determined by the s-and p-wave coalescence conditions. Our approach can be unproblematically combined with the original orbital-invariant ansatz, and can be extended to the UCCSD(F12) level of theory, using ROHF or UHF references, with minimal changes to the existing UCCSD(F12) program.

To properly assess the impact of spin-flipped geminals on the accuracy of a calculation, it is necessary to compare to analogous methods without the spin-flipped geminals. When optimising the amplitudes, the comparison is between the original orbital invariant formulation of Ref 36 and the new methods presented here. However, no open-shell fixed-amplitude method exists in the literature where the coalescence conditions are fulfilled as far as possible within the normal geminal basis, that is, without spin-flipped functions. In this work we derive such a fixed-amplitude method for MP2-F12 with UHF and ROHF references, which is both orbital invariant and size extensive and reduces to Ten-no's rational generator method when the α and β orbitals are identical.

II. THE GEMINAL BASIS FOR OPEN-SHELL F12 CALCULATIONS

Stemming from the original ideas of Kutzelnigg, 1 the geminal basis in F12 methods is chosen such that the s-and p-wave coalescence conditions could be satisfied for every spin-adapted first order pair function u s ij . The s-wave condition applies to the singlet pair functions, u 0 ij , and the p-wave condition to the triplet pair functions u 1 ij .

u 0 ij = u 0 ij (r 12 = 0) + 1 2 r 12 (|ij + |ji )(2 + 2δ ij ) -1/2 + O(r 2 12 ), (1)
u 1 ij = r 12 • ∂u 1 ij ∂r 12 r 12 =0 + 1 4 r 12 |ij + O(r 3 12 ) αα, ββ, (2)
u 1 ij = r 12 • ∂u 1 ij ∂r 12 r 12 =0 + 1 4 r 12 |ij -|ji 2 -1/2 + O(r 3 12 ) αβ. (3)
|ij is a two-electron determinant of spin-orbitals i and j, where it is understood that the first orbital has α spin and the second has β spin, for the spin-opposite case. We use this convention consistently throughout this article. Introducing variational flexibility in an orbital-invariant manner, 36 the first order pair functions are expanded as

u ij = a<b t ij ab |ab + k<l c ij kl Q12 f (r 12 )|kl , (4)
where k, l run over occupied and a, b over virtual Hartree-Fock spin-orbitals. Q12 is the strong orthogonality operator and ∂f (r 12 )/∂r 12 | r 12 =0 = 1. For restricted closed-shell systems, equations ( 1) and ( 3) can be used to select fixed F12 amplitudes.

c ij kl = 1 4 δ i k δ j l αα, ββ, (5)
c ij kl = 3 8 δ i k δ j l + 1 8 δ i l δ j k αβ. (6)
Indeed, in the limit of a large orbital basis, variationally optimised F12 amplitudes tend to these values. The derivation of equations (1-3) assumes that singlet and triplet pair functions can be constructed and thus is only rigorous for closed-shell restricted first order pair functions. For UMP2 using ROHF or UHF references, the coalescence conditions must be re-examined in terms of spin-orbitals. Our derivation of the coalescence conditions for open-shell MP2 is reported in the appendix and we consider both UHF and ROHF references.

u ij = u ij (r 12 = 0) + 1 2 r 12 (|ij + | jī )2 -1/2 +r 12 • ∂u ij ∂r 12 r 12 =0 + 1 4 r 12 |ij -| jī 2 -1/2 + O(r 2 12 ). (7)
The spin-orbital ī has the same spatial component as i, but the spin component is flipped such that ī = (ŝ + + ŝ-)i. Notably, both s-and p-wave coalescence conditions apply to the opposite spin pair functions, even though it is not possible to construct singlet and triplet functions from the corresponding pairs of occupied spin-orbitals.

u ij = a<b t ij ab |ab + k<l c ij kl Q12 f (r 12 )|kl + k< l c ij kl Q12 f (r 12 )| kl . (8)
Equation ( 7) is naturally also valid for closed shell HF references and in this case the spinflipped functions are obviously always present in the normal geminal basis and the functions f (r 12 )| kl are redundant. For ROHF references, the only spin-flipped functions missing from the normal geminal basis are those where either k or l belongs to a singly occupied restricted HF orbital and thus has no spin-flipped counterpart among the HF occupied orbitals.

A. Fixed F12 amplitudes in ROMP2-F12 and UMP2-F12

The coefficients in equation ( 8) may be selected using the coalescence conditions, which gives

c ij kl = 1 4 δ i k δ j l c ij kl =0 αα, ββ, (9)
c ij kl = 3 8 δ i k δ j l c ij kl = 1 8 δ i l δ j k αβ. ( 10
)
This method is orbital invariant and size extensive and is also size consistent when combined with a UHF reference. ). We will refer to MP2-F12 methods that use spin-flipped geminals with this choice of fixed amplitudes as MP2-F12-sp.

To assess the effect of including the spin-flipped geminals, we must define a fixed amplitude method where the geminal basis is restricted to that of equation ( 4). The best choice of fixed amplitudes is to satisfy the derivative discontinuities as far as possible within the basis available. For ROHF references and restricted orbitals, this corresponds to

c ij kl = 1 4 δ i k δ j l αα, ββ, (11)
c ij kl = 3 8 δ i k δ j l + 1 8 δ i l δ j k αβ, "d", (12)
c ij kl = 1 2 δ i k δ j l αβ, "s", (13)
where "d" denotes pairs ij where both i and j refer to orbitals that are doubly occupied and "s" pairs where either i or j is singly occupied in the ROHF reference. For "d" pairs, both the s-and p-wave discontinuities can be satisfied, but for "s" pairs, only the s-wave cusp can be satisfied. This choice differs from that of Ref. 27, where c ij kl = 3 8 δ i k δ j l was used for the "s" pairs. However, it is easily verified that in the limit of a large orbital basis, the optimised amplitudes tend to our values. Upon rotating the orbitals, for example for semi-canonicalisation or localisation, the αβ amplitudes become

c ij kl = 1 2 δ i k δ j l -1 8 mn s in s kn s mj s ml + 1 8 s il s kj , (14)
where s il is the spatial overlap between the α spin-orbital φ i and the β spin-orbital φ l .

For calculations based on UHF references, we suggest that equation ( 14) is the best choice within the constraints of the normal geminal basis, even though neither the s-nor the pwave conditions are rigorously met (the coefficients are based on the overlap and not on how closely ϕ l ϕ k matches ϕ i ϕ j at r 12 = 0). Such a choice is both orbital invariant, size extensive and, importantly, size consistent, by virtue of the dependence on the overlap s il . We will refer to this method as MP2-F12sp.

As an aside, we note that the Quantum Monte Carlo community makes a different choice for open-shell calculations, 37 which when applied to MP2-F12 theory corresponds to Here, only the s-wave cusp is satisfied for the αβ pairs. For restricted closed-shell systems, this reduces to c ij ij = 1 2 for the spin-adapted singlet pairs, and c ij ij = 1 4 , 1 2 and 1 4 for the m s = 1, 0, -1 components of the triplet excitations, respectively. Although this is both orbital invariant and size consistent, it clearly introduces spin-contamination 23 and we do not consider it further.

c ij kl = 1 4 δ i k δ j l αα, ββ, (15)
c ij kl = 1 2 δ i k δ j l αβ. (16)
Concerning the relative computational costs of the above choices, we make the following observations, which apply to density fitted MP2-F12. For closed-shell calculations, we remind the reader that the integral transformation steps scale as N 4 , that the contractions of three-index integrals to build four-index quantities scale as N 5 and that the contractions of the four-index quantities to build the B and V matrices scale as N 4 , because only the diagonal of the spin-adapted matrices is required for the fixed amplitude ansatz. For open-shell calculations:

• The N 4 steps increase by a factor of two since integrals for α and β orbitals are required. This is true for the Coulomb, exchange and Fock transformed integrals even if restricted orbitals are used. For the spin-flipped method the prefactor is increased by a further factor of two, since spin-flipped integrals are required. This further factor of two is avoided if restricted orbitals are used (RMP2-F12).

• The N 5 contractions are twice as expensive as the closed-shell case (again, this also applies to the most expensive contractions in RMP2-F12). For the spin-flipped method, a further factor of two arises due to the spin-flipped geminals, but this is avoided for RMP2-F12.

• For both the normal and spin-flipped fixed amplitude methods, N 4 scaling can be recovered for the final construction of the αβ B and V matrices by appropriately contracting four-index quantities with the fixed amplitudes. For spin-flip this is trivial and, for the overlap-based approach, due to the factorisability of the terms in equation ( 14) into separate transformations for α and β orbitals, this can be achieved by four N 5 contractions with a low prefactor (they scale at most as O 3 V 2 , where O is the number of active occupied orbitals and V is the maximum number of virtual or CABS orbitals). We note that an analogous contraction has been utilised by Köhn to simplify the equations for restricted closed-shell F12 methods. To summarise, using spin-flipped geminals is approximately twice as expensive as using the normal geminal basis with equations ( 11)-( 14), unless restricted orbitals are used, in which case the two methods have comparable cost. Numerical tests comparing the performance of the two approaches are reported in section III. ROHF-based MP2-F12, this is not necessary since the redundancies can be removed exactly by constructing the geminal basis from the restricted orbitals. However, we consider that a balanced treatment of open-and closed-shell species will only be obtained if the size of the geminal basis is proportional to the number of correlated electron pairs. We therefore propose that equation ( 8) should be replaced by

u ij = a<b t ij ab |ab + k<l c ij kl Q12 f (r 12 )(c|kl + c| lk ), (17)
where the contraction coefficients c and c are determined using the coalescence conditions

c = 1 4 c =0 αα, ββ, (18)
c = 3 8 c = 1 8 αβ. (19)
It should be noted that for same-spin pairs, the spin-flipped geminals always have exactly zero contribution. Variational optimisation of the c ij kl is no more complicated than if equation (4) is used. Resticting c ij kl to the diagonal amplitudes c ij ij results in a non-orbital invariant method and we do not consider it here, although it may be useful in combination with localisation techniques. In section III we present results of calculations where the F12 amplitudes are optimised with and without the flipped contributions. We will refer to the standard orbital invariant method, without spin-flipped geminals, as MP2-F12-õ and to the new method The relative computational cost of including or excluding spin-flipped geminals has been discussed in the previous section for the fixed-amplitude MP2-F12 methods. If the coefficients in equation ( 4) or (17) are variationally optimised, then N 6 contractions are required for the construction of the B and V matrices, and also for the inclusion of coupling terms between t ij ab and c ij kl when solving the amplitude equations. Extending the approach of Bokhan et al. to the invariant ansatz for UHF-based MP2-F12 would lead to a fourfold increase in the N 6 contractions for the αβ B matrix compared to neglecting spin-flipped geminals and a twofold increase for the αβ V matrix. If an ROHF reference is used, this corresponds to the MP2-F12-o + method and only the spin-flipped geminals for the "s" pairs need be computed in additional to the the normal geminals, leading to a slight increase in the prefactor for the N 6 contractions compared to MP2-F12-õ. However, if the geminal basis is contracted as we suggest, there is no additional computational cost for including the spin-flipped functions for the N 6 steps.

C. The geminal basis in open-shell CCSD-F12 calculations

At the CCSD level of theory it has been standard practice 38 to use the geminal basis appropriate for MP2 and also to use the MP2 coalescence conditions to define fixed F12 amplitudes, 28,29,32 even though the MP2 and coupled-cluster coalescence conditions differ. 39 Köhn has recently included additional geminals that amend this deficiency, demonstrating that they are important for response properties. 33 For ground state energies, however, selecting the geminal basis and amplitudes based on MP2 appears to be sufficient and we pursue this approach here. Switching to second quantisation, the doubles cluster operator becomes

T2 = 1 4 aibj t ij ab τ ab ij + 1 8 kilj c ij kl µν w kl µν τ µν ij . ( 20
)
where µ, and ν are orbitals in the complete basis and

w µν kl = c µν| Q12 f (r 12 )|kl + c µν| Q12 f (r 12 )| lk (21)
The contraction coefficients are c = 11) and ( 14) are recommended in the absence of spinflipped geminals and if spin-flipped geminals are included, then c ij kl = δ i k δ j l . If a UCCSD-F12 program is available that uses the normal geminal basis, all that is required to extend it to employ the geminal basis is to replace the F12 integrals with the above linear combination of normal and flipped integrals. In section IV we present results of spin-flipped calculations at the ROHF-based UCCSD(F12) level of theory.

III. MP2-F12

Our focus is the effect of the spin-flipped geminals on the correlation treatment. For ROHF-based MP2, contributions from both single and double excitations occur, but they are uncoupled and the singles energy depends only on the Fock matrix elements. Since the basis set convergence for the singles and doubles are thus very different, we partition the MP2 energy into HF + = E HF + f a i t i a and the doubles correlation contribution.

A. Atomisation energies

We have recently computed a set of basis set limit frozen-core (fc) MP2 energies for a set of 106 molecules of H, C, N, O and F to an accuracy of 99.95%, which are intended as a test set for F12 methods. 40 Here we use these benchmarks to assess the relative performance of the MP2-F12sp, MP2-F12-sp, MP2-F12-õ, MP2-F12-o and MP2-F12-o + methods for atomisation energies using the cc-pVXZ-F12 orbital basis sets 20 with X = D, T and Q. All UMP2-F12 calculations are based on RHF or ROHF references and were performed using the Turbomole program package using the default options, namely Ansatz 2, approximation B, the specially optimised complementary auxiliary basis sets of Peterson and co-workers 21 and the aug-cc-pwCV(X+1)Z MP2-fitting 41 and aug-cc-pV(X+1)Z RI-JK basis sets 42 for density fitting. The Slater-type correlation factor was expanded in six Gaussians 43 and the exponents for the Slater function were 0.9, 1.0 and 1.1 a -1 0 for the cc-pVXZ-F12 basis sets, X = D, T and Q, and 1.4 a -1 0 for the aug-cc-pVXZ basis sets. The integrals involving the single commutator of f (r 12 ) with the kinetic energy were computed from the matrix representation of the Fock and exchange matrices. Furthermore, a singles perturbative correction was added to the HF energy using equation (48) of Ref. 27 shell case). In our previous work, the basis set limits for the atoms were computed using the MP2-F12-õ method. The corresponding MP2-F12-o energies are slightly closer to the basis set limit and we use these improved benchmarks (Table I) in this work.

In Table II we report the mean ∆ and standard deviation σ (s N ) for basis set errors of atomisation energies using MP2-F12 methods with and without the spin-flipped geminal basis. Results for both the fixed amplitude (MP2-F12sp, MP2-F12-sp) and optimised (MP2-F12-õ, MP2-F12-o, MP2-F12-o + ) approaches are presented. In all cases, the inclusion of spin-flipped geminals decreases both the mean and standard deviation of the basis set errors for the correlation energy contributions, which indicates that a more balanced treatment of atoms and molecules has been achieved. For the triple-and quadruple-zeta quality basis sets, the improvement is at least a factor of two.

Comparing the MP2-F12sp and MP2-F12-sp results, we note that if the spin-flipped geminals are excluded, the correlating basis set for the atoms is poorer than that for the molecules, contrary to all experience with conventional MP2 calculations. As a result, the basis set errors for the HF and correlation contributions have on average opposite signs, which leads to a cancellation of errors for the total MP2 atomisation energies. By including the spin-flipped geminals, the correlating basis set for the atoms is better than that of the molecules, as usually expected, which reverses the sign of the errors and results in a deterioration in the mean basis set errors for the total MP2 atomisation energies (not reported in Table II). The improvement in the standard deviation, however, is retained.

For the methods where the c ij kl are variationally optimised, the mean basis set errors are smaller in magnitude than for the fixed-amplitude methods, while the standard deviations are similar. Since the errors for the spin-flipped and normal geminal calculations decrease from opposite directions, this represents an improvement in the correlation description rather than a cancellation of errors due to basis set superposition error, which would favour the molecules over the atoms. The effect of fully optimising both the normal and spin-flipped geminal amplitudes may be seen by comparing the results for MP2-F12-o and MP2-F12o + . We find that the basis set errors are fairly similar and conclude that the contraction coefficients c and c used in the MP2-F12-o method are close to optimal. The mean error for MP2-F12-o + is more negative than that of MP2-F12-o, which indicates that the geminal basis is less balanced, favouring the atoms over the molecules. We therefore recommend the MP2-F12-o approach.

B. Ionisation potentials and electron affinities

The statistics for the atomisation energies in the previous section give an indication of how balanced the method is for treating open and closed shell systems. Here we investigate the effect of spin-flip on the relative accuracy of atoms and ions. In Table III we report basis set errors of ionisation potentials and electron affinities of C, N, O and F, computed using fc-MP2-F12/cc-pVXZ, X = D, T and Q, with and without spin-flipped functions. The effect of spin-flip on the basis set error of the total correlation energy is largest for the s 2 p 4 configuration, and the importance of spin-flip decreases monotonically when progressing to either the s 2 p 0 or the s 2 p 6 configurations. Thus the IPs for C, N and O are improved with spin-flip, while F deteriorates. Similarly, the EAs of O and F deteriorate. The mean absolute basis set error over the seven energy differences does not improve upon incorporating spinflipped geminals into the wave function, but the errors do become more uniform.

Clearly, when computing electron affinities, the capability of the basis to represent the diffuse electron cloud is an issue. The cc-pVXZ-F12 basis sets do not appear to contain sufficiently diffuse functions for the correlation treatment of the anions. Moreover, our investigations show that the aug-cc-pVXZ basis sets do not perform better than the cc-pVXZ-F12 basis sets for the EAs computed here. It should be noted that the HF error has been greatly reduced by the CABS singles correction. The necessity for more diffuse functions is particularly important in F12 methods, which may be seen from the dependence of the basis set errors on the exponent γ of the correlation factor. Figures 123show the γ dependence of the MP2-F12-sp and MP2-F12-o basis set errors of C, O and F, respectively, together with their cations and anions, using the cc-pVTZ-F12 basis. The MP2-F12-o curves for the anions differ strongly from those of the atoms and cations and lie far below the corresponding MP2-F12-sp curves (for O and F) at low values of γ. The exaggerated importance of off-diagonal (and diagonal for C) excitations into long-range geminals for the MP2-F12-o method is due to the lack of sufficiently diffuse functions in the orbital basis and leads to the undesirable γ dependence. Figures 123also demonstrate that choosing γ based on the neutral atom necessarily biases the geminal basis in favour of the atom. The exponent γ is related to the extent of the correlation hole, which is naturally smaller for the cation (larger γ) with a contracted electron density and larger for the anion (smaller γ). In the case of the cations, however, this bias does not appear to be severe. Benchmark fc-CCSD correlation energies for a subset of 30 of the 106 molecules of our test set have recently been computed using the CCSD(F12) method with the same orbital and auxiliary basis sets as for the fc-MP2 benchmarks. 45 The corresponding benchmark valence CCSD correlation energies for the atoms are reported in Table I, computed using the UCCSD(F12)-o method with semi-canonical ROHF orbitals. We expect that just as for the MP2 correlation energies, these benchmarks are accurate to 99.95% and we used them to assess the effect of including spin-flipped functions in CCSD(F12) calculations. In Table IV we present basis set errors for the CCSD(F12)sp, CCSD(F12)-sp, CCSD(F12)-õ and CCSD(F12)-o methods for atomisation energies using the cc-pVXZ-F12 orbital basis sets. The same auxiliary basis sets and F12 ansatz and approximations were used as for the MP2-F12 calculations discussed in section III A.

We find that the improvements obtained for MP2-F12 by including spin-flipped geminals are mirrored in the CCSD(F12) results and the discussion for MP2-F12 applies equally well to CCSD(F12). We note, however, that the standard deviation of the basis set error is somewhat larger on average for CCSD than for MP2. One possible source for this additional uncertainty is the incomplete basis for the singles contribution to the correlation energies, which is more important for open-shell systems than closed-shell systems. Unfortunately, the decomposition of the basis set limit fc-CCSD energies into doubles and singles contributions is not available for the molecules and we are unable to examine these separately. This is possible for our current calculations on atoms and ions and we investigate the singles and doubles contributions separately for IPs and EAs in the next section.

B. Ionisation potentials and electron affinities

In Table V we present the basis set errors for ionisation potentials and electron affinities for the CCSD(F12)sp, CCSD(F12)-sp, CCSD(F12)-õ and CCSD(F12)-o methods, computed using fc-UCCSD(F12) with semi-canonical ROHF orbitals. We partition the CCSD(F12) energy into HF + and the doubles correlation contribution, where HF + is the sum of the HF energy, the singles correlation, and the CABS singles correction. Just as for the atomisation energies, we find that the effect of spin-flip on the CCSD energies is very similar to that of the MP2 energies, discussed in section III B. Comparing the HF + and doubles components of the CCSD energies with the MP2 equivalents, we find that the singles contributions converge at similar rates. The basis set errors from the doubles are larger on average for CCSD than MP2 for the cc-pVDZ-F12 basis, but are similar to MP2 for the cc-pVQZ-F12 basis. This is consistent with the effect of the singles basis set incompleteness on the doubles, through the coupling which is present in CCSD but not MP2.

V. CONCLUSION

For closed-shell calculations, the geminal basis functions are chosen as f (r 12 )|ij , where i j are spin-orbitals, occupied in the HF reference. These are sufficient to exactly satisfy the s and p wave MP2 coalescence conditions and thus accelerate the basis set convergence from X -3 to X -7 . As first noted by Bokhan et al., spin-flipped geminal functions are required to achieve this for open-shell calculations. In this work, we have defined a hierarchy of open-shell F12 methods with and without spin-flipped functions which we have used to investigate the numerical importance of the spin-flipped functions for atomisation energies, ionisation potentials and electron affinities at the MP2 and CCSD levels of theory. In particular, we have solved the issue of redundancy or near linear dependency in the normal and spin-flipped basis by using contracted geminals with contraction coefficients determined from the coalescence conditions. Furthermore, we have used an analysis of the spin-orbital MP2 coalescence conditions to establish the correct choice of fixed amplitudes when no spin-flipped functions are used in the calculation. All of the methods we present are orbital invariant, size extensive and size consistent when combined with UHF references and reduce to either the sp or original orbital-invariant ansatz for the closed-shell case.

Our calculations on a set of 106 molecules of H, C, N, O and F convincingly demonstrate that spin-flipped functions are important for atomisation energies, where a balanced treatment of open-and closed shell-species is necessary. The basis set errors are approximately halved by including the spin-flipped functions for both MP2 and CCSD calculations. Furthermore, we find that our contracted geminal approach is superior to optimising the amplitudes for both the normal and spin-flipped geminals, since the number of variational parameters is then proportional to the number of correlated electron pairs. For ionisation potentials and electron affinities, we find that including spin-flipped geminals improves the consistency of the energies, but that the accuracy is limited due to difficulties in selecting an appropriate exponent for the correlation factor and the need for diffuse functions for electron affinities.

In conclusion, we find that spin-flipped geminals are very important for open-shell F12 calculations and recommend that they be used. For the fixed amplitude approach, the sp ansatz should be used as proposed by Bokhan et al. in Ref. 35 and if the F12 amplitudes are to be optimised, our contracted geminal approach should be used. If the MP2 method is desired and computational cost is an issue, then our sp method is recommended. For coupled-cluster calculations, the additional cost of incorporating them in either the sp or optimised amplitude methods is only at the level of the F12 integrals and is thus minimal.

VI. APPENDIX: UMP2 COALESCENCE CONDITIONS

The exact solution to the Born-Oppenheimer, non-relativistic Schrödinger equation may always be expressed as a sum of products of spatial and spin eigenfunctions. Correspondingly, Kato's proof of the cusp condition at electron coalescence 24 and the subsequent analysis of Pack and Byers Brown, 25 were concerned with spatial eigenfunctions only. For approximate wave functions that are not eigenfunctions of spin, such as UMP2, the analysis leading to the coalescence conditions must be repeated in terms of spin-orbitals. Fortunately, this is relatively straightforward, at least at the level of UMP2, and follows the analysis in Refs.

32 and 39 for RHF references rather closely.

Consider canonical or semi-canonical orbitals, which are exact eigenfunctions of the zeroth-order operator. The equation for the first-order wave function separates into equations for each pair,

V1 V2 f 0 1 + f 0 2 -ε i -ε j |u ij + V1 V2 1 r 12 |ij = 0. ( 22
)
15 Here

|u ij = (1 -P12 )u ij / √ 2 and |ij = (1 -P12 )φ i φ j / √ 2
, where φ i and φ j are occupied spin-orbitals and P12 is the permutation operator. V1 V2 is the strong orthogonality projector

(1 -Ô1 )(1 -Ô2 ), where Ô projects onto the space of occupied spin-orbtals. For the ROHF case, f 0 = ff 1 , where f 1 = Ô f V + V f Ô and the f 1 terms drop out since they do not connect the space V 1 V2 with |u ij . For the UHF case, f 0 is the usual Fock operator ( f 1 = 0).

The presence of f 1 for ROHF leads to an equation for single excitations in addition to equation (22), but this does not effect the nature of the first-order wave function at electron coalescence and we do not consider it further.

The analysis of Pack and Byers Brown proceeds by transforming to centre of mass, s = (r 1 + r 2 )/2, and relative coordinates, r 12 = r 2r 1 and equating powers of r 12 as far as possible. To manipulate equation ( 22) into a form amenable to such an analysis, we follow Ref. 39 and introduce a pair function w ij that has no orthogonality constraints, such that

u ij = V1 V2 w ij .
For UHF, [ f , V ] = 0 and we may drop the strong orthogonality projector from the equation. For ROHF [ f , V ] = Ô f V -V f Ô. Finally, we remove spin by left projecting onto the spin function αα, ββ or αβ. In place of equation ( 22), we now have

f1 + f2 -ε i -ε j ω ij + 1 r 12 Θ ij = 0, (23)
where f is replaced by ff Ô for ROHF. Note that for same-spin pairs, P12 Θ ij = -Θ ij and P12 ω ij = -ω ij , but for spin-opposite pairs, Θ ij and ω ij are just the spatial components of φ i φ j and w ij , respectively. Expanding the Fock operator f in powers of r 12 , we obtain

f1 + f2 = -∇ 2 r 12 + O(r 0 12 )Y 00 + O(r 1 12 ), (24)
where Y lm are the surface harmonics for the angular dependence of r 12 in spherical coordinates. This differs from the closed shell case in that the O(r 1 12 ) terms, which are angle dependent, do not cancel. To obtain this, we have used

1 r i1 = 1 r is + r 12 r 2 is cos θ i + O(r 2 12 ), ( 25
) 1 r i2 = 1 r is - r 12 r 2 is cos θ i + O(r 2 12 ), ( 26
)
where r is is the distance between s and the centre r i and cos θ i is the cosine between the vectors r is and r 12 . For ROHF, the presence of f Ô does not change equation (24). The spatial components of the zeroth-and first-order pair functions are now expanded in radial and

where the second term is a spin-flipped determinant, for example, with the alpha spin-orbital φj = ϕ j α, where ϕ j is the spatial component of the beta spin-orbital φ j . From equation (32) it is clear that any linear combination of ϕ i (r 1 )ϕ j (r 2 ) and ϕ j (r 1 )ϕ i (r 2 ) can be used to satisfy the s-wave cusp and equations ( 35) and ( 36) are not unique solutions. However, the above choice for equation (36)

  f (r 12 ) is the correlation factor. To satisfy the s-and p-wave coalescence conditions for open-shell MP2-F12 calculations, Bokhan et al. introduced a spin-flipped geminal basis and

  In complete agreement with the work of Bokhan et al., we find that equation (2) holds for the same-spin pairs and for the αβ pairs

33

 33

  B. Optimised F12 amplitudes in ROMP2-F12 and UMP2-F12 Equation (4) presents no special difficulties upon variational optimisation of the F12 amplitudes c ij kl , beyond those sometimes encountered for closed-shell calculations. Indeed, this method has been employed for open-shell R12 calculations for the past 15 years. Equation (8), on the other hand contains large redundancies between the normal and spin-flipped geminals for restricted open-or closed-shell HF references, which become severe near linear dependencies when an unrestricted HF reference is used. Bokhan et al. have proposed to delete geminals based on singular value decomposition to remove the linear dependence. For

8

 8

  the flipped functions as MP2-F12-o. For comparison, we also report ROHF-based calculations where the non-redundant set of normal and flipped variational parameters are optimised separately, which will be denoted MP2-F12-o + .

1 4 ,

 4 c = 0 for like-spin pairs and c = 3 8 , c =1 8 for αβ pairs, or simply c = 1 and c = 0 if spin-flipped geminals are not to be included. For the method, equations (

13

 13

FIG. 1 FIG. 2 FIG. 3 FIG. 3 :

 1233 FIG.1The γ dependence of MP2-F12-o (dashed) and MP2-F12-sp (solid) basis set errors (meV) in the doubles correlation energies for C + (×), C (•) and C -(⋄), using the cc-pVTZ-F12 orbital basis.

  for ROHF (see Ref. 30 or 44 for the closed-

  Both the s-and p-wave coalescence conditions are met for the αβ pairs if

				Molecular Physics	Page 16 of 28
	angular functions of r 12 , with a Taylor expansion for the radial functions, whose coefficients
	depend on s.	|u ij =	3 8	V1 V2 r 12 |ij +	1 8	V1 V2 r 12 | jī .
		Θ ij (r 12 , s) =
	F o r F o r
		P P
		e e
		e r e r
				R R
				e e
				v i e v i e
					w w
						O n l O n l
						y y

URL: http://mc.manuscriptcentral.com/tandf/tmph

  corresponds to the correct formula for the case of identical α and β spatial orbitals, where spin and space symmetry are properly coupled. The choice in equation(35) does not require any spin-flipped functions.

	F o r F o r
	P P
	e e
	e r e r
	R R
	e e
	v i e v i e
	w w
	O n l O n l
	y y

TABLE I :

 I Basis set limit HF and valence MP2 and CCSD correlation energies for the atoms C, N, O and F in E h .

	Atom	E HF	MP2 correlation CCSD correlation
	C carbon	-37.6886122	-0.0818650	-0.0980777
	N nitrogen	-54.4009236	-0.1138167	-0.1258611
	O oxygen	-74.8093817	-0.1792901	-0.1910343
	F fluorine F o r	-99.4093241	-0.2477734	-0.2534840
	P
	e
		e r
		R
		e
			v i e
			w
			O n l
				y

TABLE II :

 II Basis set error statistics for the HF + and doubles correlation contributions to fc-MP2-F12 atomisation energies (kJ/mol per valence electron). HF + refers to E HF + f a i t i a + f p ′

				i t i p ′ .
		cc-pVDZ-F12 cc-pVTZ-F12 cc-pVQZ-F12
	MP2-F12	∆ σ	∆ σ	∆ σ
	HF + sp	-0.13 0.05 0.16 0.23	-0.02 0.01 0.08 0.06	0.00 0.00 0.03 0.02
	sp õ o o + F o r	-0.34 0.11 0.09 0.12 -0.14 0.08 -0.18 0.09	-0.07 0.03 0.07 0.04 -0.02 0.02 -0.05 0.03	-0.02 0.01 0.03 0.02 -0.01 0.01 -0.01 0.01
	P
	e
		e r
		R
		e
		v i e
			w
			O n l
				y

TABLE III :

 III Basis set errors in the fc-MP2-F12 HF + and doubles correlation contributions to ionisation potentials and electron affinities (meV). HF + refers to E HF + f a i t i a + f p ′ i t i p ′ .

				IP				EA
	Basis	MP2-F12 C	N	O	F	C	O	F	MAD
	cc-pVDZ-F12 HF +	5.9	6.4	3.5	4.2	1.7	1.4	1.1	3.5
		sp sp	10.6 8.4	10.8 8.3	59.0 40.1	2.3 39.4	40.1 40.8	24.1 52.1	-22.5 45.7	24.2 33.5
	õ o cc-pVTZ-F12 HF + sp sp F o r 11.1 8.1 1.6 4.1 3.0	12.5 9.2 2.0 3.8 2.7	32.3 28.8 2.3 15.1 10.6	12.9 24.5 2.9 -1.8 7.7	39.8 38.5 0.7 28.9 28.8	43.0 52.8 1.4 21.6 29.9	26.9 46.4 1.7 4.8 21.1	25.5 29.8 1.8 11.4 14.8
		õ o	4.3 3.2 P	4.3 3.1	11.1 8.9	2.0 6.2	29.1 28.5	26.4 30.5	13.9 22.4	13.0 14.7
	cc-pVQZ-F12 HF + sp sp õ	0.2 e 0.2 1.1 1.3 0.7 1.0 1.2 1.4 e r	0.5 4.8 3.3 4.1	0.4 -0.4 2.2 0.3	0.0 12.9 12.9 12.9	0.0 22.5 24.7 23.4	-0.4 9.4 14.1 11.0	0.2 7.5 8.4 7.8
		o	0.7	1.1 R 3.2	2.0	12.7	24.8	14.3	8.4
				e
					v i e
						w
						O n l
								y

TABLE IV :

 IV Basis set error statistics for the fc-CCSD(F12) correlation contribution to atomisation energies (kJ/mol per valence electron).

	cc-pVDZ-F12 cc-pVTZ-F12 cc-pVQZ-F12
	CCSD(F12)	∆ σ	∆ σ	∆ σ
	sp	0.37 0.59	0.12 0.14	0.07 0.07
	sp õ	-0.24 0.33 0.23 0.43	-0.06 0.07 0.08 0.07	0.02 0.04 0.09 0.04
	o	0.04 0.37	0.03 0.04	0.07 0.03
	F o r
	P
	e
	e r
		R
		e
		v i e
			w
			O n l
				y

TABLE V :

 V Basis set errors in the CCSD(F12) HF + and correlation contributions to ionisation potentials and electron affinities (meV).

				IP				EA
	Basis	CCSD(F12) C	N	O	F	C	O	F	MAD
	cc-pVDZ-F12	HF + sp sp	6.5 13.3 7.1	6.5 20.6 11.5	11.7 83.5 69.2	4.8 16.5 59.4	0.9 42.3 38.6	1.7 30.3 64.0	-8.2 -22.1 51.6	5.7 32.6 43.1
		õ	16.3	24.2	51.5	22.8	48.6	44.0	20.3	32.5
		o	12.9	17.4	49.6	34.5	44.7	54.2	37.6	35.8
	cc-pVTZ-F12	HF +	1.7	2.1	2.6	3.2	0.4	3.1	0.6	2.0
		sp	3.4	4.5	23.0	4.3	28.7	23.1	4.4	13.1
		sp	1.5	2.1	19.0	15.0	27.7	33.0	22.4	17.2
		õ	3.6	4.7	12.4	6.5	30.0	25.7	13.4	13.8
		o	2.9	3.4	10.6	9.4	29.4	28.7	19.8	14.9
	cc-pVQZ-F12	HF + sp	0.3 0.8	0.3 1.5	0.9 8.8	0.5 1.9	-0.2 12.8	0.9 22.0	-1.0 8.1	0.6 8.0
		sp	0.2	0.7	7.4	5.0	12.5	24.6	13.6	9.1
		õ	0.6	1.3	3.1	0.5	13.0	20.4	8.6	6.8
		o	0.4	0.9	2.6	1.6	12.8	21.3	10.7	7.2

Acknowledgments

This research has been supported by the Deutsche Forschungsgemeinschaft through grant nos. TE 644/1-1 and KL 721/2-3.

Note that for same-spin pairs, only terms with odd l survive, but for αβ pairs, all l are present. Following Pack and Byers Brown, equation (23) becomes

Equating terms in r -1 12 gives the familiar s-and p-wave coalescence conditions (l = 0 and 1, respectively),

For the same-spin pairs, the p-wave conditions apply. This result is of course well known.

For spin-opposite pairs, both the s-and p-wave conditions apply so that

Since the Y 00 terms are symmetric with respect to the operator P12 , we may identify

where ϕ i is the spatial component of φ i . Similarly, since the Y 1m terms are antisymmetric, we may identify

For the same-spin (triplet) pairs, the p-wave coalescence conditions are satisfied if

For the αβ pairs, the s-wave cusp conditions are satisfied if


Ans M W Van Den Ouweland
email: a.vandenouweland@erasmusmc.nl

Peter Elfferich

Bernard Zonnenberg

Willem F Arts

Tjitske Kleefstra

Mark Nellist

Jose M Millan

Caroline Withagen- Hermans

Anneke J A Maat -Kievit

Dicky J J Halley

Characterisation of TSC1 promoter deletions in Tuberous Sclerosis Complex patients

Keywords: TSC1, promoter, deletions, MLPA

Tuberous Sclerosis Complex (TSC), an autosomal dominant disorder, is a multisystem disease with manifestations in the central nervous system, kidneys, skin and/or heart. Most TSC patients carry a pathogenic mutation in either TSC1 or TSC2.

All types of mutations, including large rearrangements, nonsense, missense and frameshift mutations, have been identified in both genes, although large rearrangements in TSC1 are scarce.

Here we describe the identification and characterisation of 8 large rearrangements in TSC1 using Multiplex Ligation-dependent Probe Amplification (MLPA) in a cohort of 327 patients in whom no pathogenic mutation was identified after sequence analysis of both TSC1 and TSC2 and MLPA analysis of TSC2. In 4 families, deletions only affecting the noncoding exon 1 were identified. In one case, loss of TSC1 mRNA expression from the affected allele indicated that exon 1 deletions are inactivating mutations. Although the number of TSC patients with large rearrangements of TSC1 is small, these patients tend to have a somewhat milder phenotype compared to the group of patients with small TSC1 mutations.

INTRODUCTION

Tuberous Sclerosis Complex (TSC, MIM#191100) is an autosomal dominant disorder characterised by seizures, mental retardation and hamartomas in multiple organ systems, including brain, skin, heart, lungs and kidneys [START_REF] Crino | The tuberous sclerosis complex[END_REF] . Mutations in either TSC1 or TSC2 are the underlying cause of the clinical symptoms in TSC patients. In about 75 -85% of the patients meeting the definite clinical criteria a pathogenic TSC1 or TSC2 mutation is identified [START_REF] Au | Genotype/phenotype correlation in 325 individuals referred for a diagnosis of tuberous sclerosis complex in the United States[END_REF][START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF][START_REF] Dabora | Mutational analysis in a cohort of 224 tuberous sclerosis patients indicates increased severity of TSC2, compared with TSC1, disease in multiple organs[END_REF][START_REF] Jones | Application and evaluation of denaturing HPLC for molecular genetic analysis in tuberous sclerosis[END_REF][6][START_REF] Van Slegtenhorst | Identification of the tuberous sclerosis gene TSC1 on chromosome 9q34[END_REF] . The genes are categorised as tumour suppressor genes, since loss of heterozygosity has been shown in TSC-associated lesions [START_REF] Henske | Allelic loss is frequent in tuberous sclerosis kidney lesions but rare in brain lesions[END_REF] .

TSC1 consists of 23 exons, of which exon 1 and 2 are noncoding. A core promoter has been defined by functional analysis [START_REF] Ali | Identification of a core promoter and a novel isoform of the human TSC1 gene transcript and structural comparison with mouse homolog[END_REF] . This region of 587 bp of size is situated 510 bp upstream of exon 1 and runs into exon 1. No TATA or CAAT boxes are present in this promoter region. Several transcription factor binding sites are present including SP1, E2F and GATA sites. For the detection of small (point) mutations in TSC1 and TSC2, several screening technologies have been undertaken: Denaturing Gradient Gel Electrophoresis (DGGE), Single Strand Conformation Polymorphism (SSCP), Protein Truncation Test (PTT), Denaturing High Pressure Liquid Chromatography (DHPLC) and direct sequencing [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF][START_REF] Jones | Application and evaluation of denaturing HPLC for molecular genetic analysis in tuberous sclerosis[END_REF][START_REF] Mayer | Mutation screening of the entire coding regions of the TSC1 and the TSC2 gene with the protein truncation test (PTT) identifies frequent splicing defects[END_REF][START_REF] Van Bakel | Mutations in the TSC2 gene: analysis of the complete coding sequence using the protein truncation test (PTT)[END_REF][START_REF] Rendtorff | Analysis of 65 tuberous sclerosis complex (TSC) patients by TSC2 DGGE, TSC1/TSC2 MLPA, and TSC1 long-range PCR sequencing, and report of 28 novel mutations[END_REF][13][START_REF] Roberts | Denaturing highperformance liquid chromatography (DHPLC) is a highly sensitive, semiautomated method for identifying mutations in the TSC1 gene[END_REF] . Larger rearrangements have been detected by Fluorescence in situ hybridisation (FISH), Southern blotting, Long Range (LR)-PCR and Multiplex Ligation-dependent Probe Amplification (MLPA) analysis [START_REF] Peral | Deletion of the TSC2 and PKD1 genes associated with severe infantile polycystic kidney disease--a contiguous gene syndrome[END_REF][START_REF] Verhoef | Somatic mosaicism and clinical variation in tuberous sclerosis complex[END_REF][START_REF] Dabora | Characterisation of six large deletions in TSC2 identified using long range PCR suggests diverse mechanisms including Alu mediated recombination[END_REF][START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] . Mutations in TSC2 are more common than in TSC1, particularly in sporadic cases. Interestingly, whereas large rearrangements account for approximately 10% of all TSC2 mutations identified to date, they appear to be much less frequent in TSC1. To our knowledge, only 8 different TSC1 deletions have been described so far [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF][START_REF] Longa | TSC1 and TSC2 deletions differ in size, preference for recombinatorial sequences, and location within the gene[END_REF][START_REF] Nellist | Large deletion at the TSC1 locus in a family with tuberous sclerosis complex[END_REF] . MLPA analysis of TSC1 was undertaken in patients suspected of TSC in whom no pathogenic mutation had been identified in either TSC1 or TSC2. In 4 cases a deletion of the noncoding exon 1 was identified and in a further 4 cases multi-exon deletions were detected. The deletions were characterised and it was demonstrated that deletion of exon 1 prevents TSC1 expression.

MATERIAL AND METHODS

Patient samples

Samples of patients with either a putative or definite clinical diagnosis of TSC were received for mutation analysis. Details on clinical symptoms were obtained from the referring physician using a standardised clinical evaluation form [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] .

Mutation analysis

Extraction of DNA from peripheral blood cells was performed according to standard techniques. Mutation analysis of TSC1 and TSC2 was performed by DGGE [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] or by direct sequence analysis of all coding exons and exon/intron boundaries (primers available on request). For the detection of large rearrangements in TSC2, Southern blotting, FISH and/or MLPA were performed. After the introduction of MLPA for TSC1, all patients without an identified pathogenic mutation were tested using the SALSA MLPA kit P124 (MRC-Holland, Amsterdam, The Netherlands).

MLPA was performed according to the manufacturer's instructions; products were run on an automated sequencer (ABI 3730XL, Applied Biosystems, Foster City, CA, USA) and data were analysed using Genemarker version 1.5 (Softgenetics, State College, PA, USA). If possible, all pathogenic mutations were confirmed on an independent DNA sample.

Quantitative (Q)-PCR, LR-PCR and sequence analysis of breakpoints.

All apparent deletions detected by MLPA were confirmed by further delineation of the breakpoint regions using Q-PCR, followed by LR-PCR and sequence analysis.

Real-time Q-PCR was performed using Fam-labelled Taqman assays [START_REF] Bruge | A novel Real Time PCR strategy to detect SOD3 SNP using LNA probes[END_REF] .

Primers were designed with Primer Express 2.0.0 (Applied Biosystems) in the vicinity of single nucleotide polymorphisms (SNPs) mapping to the TSC1 locus (Table 1).

Primer specificity was checked by performing BLAST analysis. Taqman probes were synthesised with a melting temperature (Tm) 8-10˚C higher than the primers by incorporating Locked Nucleic Acid (LNA) monomers in the probe. Tm values for the LNA probes were calculated using the Exiqon website (http://lna-tm.com/). The LNAbased Taqman assays were manufactured by Eurogentec (Maastricht, The Netherlands).

Gene dosage alterations were detected on an ABI7500 Real time PCR system (Applied Biosystems) by performing a relative quantification run. Real time PCR reactions were performed in a total volume of 25 μl, containing 20 ng DNA, 1 x qPCR mastermix Plus -low ROX (Eurogentec: RT-QP2x-03-WOULR), 1 x RNAse P (endogenous control) (Applied Biosystems), 30 μM forward and reverse primers and 10 μM probe. PCR conditions were as follows: an initial 2 min incubation at 50˚C, followed by 95˚C for 10 min and then 40 cycles of 95˚C for 15 sec and 60˚C for 1 min. All samples were analysed in triplicate and compared to a normal control sample LR-PCR was performed with the Expand Long Template PCR System (Roche Applied Science, Indianapolis, IN, USA). LR-PCR products were sequenced using an automated sequencer (ABI 3730XL). Nomenclature of the deletions is according to the recommendations of the Human Genome Variation Society, using reference sequence NM_000368 (17 Dec 2004; build 36, NCBI).

RNA analysis

Fibroblasts were cultured according to standard procedures. To increase the probability of recovering (truncating) mutant TSC1 RNA, nonsense mediated decay of RNA was prevented by adding cycloheximide to the cells 4.5 hours before harvesting.

RNA was isolated using the RNeasy Mini kit (Qiagen Inc., Valencia, CA, USA).

Reverse transcriptase (RT)-PCR (oligo-dT primed) was performed using the Omniscript ReverseTranscription kit (Qiagen). The primers used for RNA analysis were: Exon 20 Forward: 5'-TGTAAAACGACGGCCAGTACAGGCAGCTGTTGGTTCTT-3' Exon 23 Reverse: 5'-CAGGAAACAGCTATGACCGCCAGATGCCTCTTCATTGT-3' Exon 20/21 Forward: 5'-TGTAAAACGACGGCCAGTGCACTCAGATACCACAAAGGAA-3' Exon 23 Reverse: 5'-CAGGAAACAGCTATGACCTCTGAGCACCCGTCATTACA-3' A first round PCR was performed, followed by a nested PCR using 1μl of the first round PCR product. The PCR conditions were: 10 sec at 94˚C, followed by 10 cycles of 30 sec 94˚C, 30 sec 68˚C, 1 min 72˚C with a decrease of 1˚C in the annealing temperature per cycle, an additional 25 cycles with an annealing temperature of 58˚C and finally 5 min at 72˚C and 5 min at 20˚C. PCR products were directly sequenced using an automated sequencer (ABI 3730XL). Data were analysed using the SeqScape software (version 2.6; Applied Biosystems).

RESULTS

Mutation analysis

Mutation analysis was performed in a diagnostic cohort of 986 TSC cases. This cohort includes the group (n= 362) previously described [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] . Of those 362 patients, 276 (76%) had a definitive clinical diagnosis of TSC. In the total group of 986 patients it was not possible to give a percentage of patients meeting the clinical diagnosis TSC, since no clinical information was available of the new patients. In 172 cases (17.4%), a pathogenic mutation in TSC1 was identified, whereas TSC2 mutations were present in 487 cases (49.3%) (data not shown). In 327 cases (33.2%) no pathogenic mutation was identified in TSC1 (by direct sequence or DGGE analysis of all coding exons) or TSC2 (by direct sequence, DGGE, Southern, FISH and MLPA analysis). MLPA analysis of TSC1 in these 327 patients showed abnormal patterns in 8 unrelated patients: in 4 cases (patient numbers 30628, 21722, 21899 and 1264; Figure 1B -E) a deletion of the noncoding exon 1 was detected, 1 patient (31457; Figure 1F) had a deletion of exons 2 -23, 2 patients (29445 and 28121; Figure 1G -H) had a deletion of exons 9 -23 and one patient (14249; Figure 1I) was identified with a total gene deletion (Figure 1).

Characterisation of the breakpoints

Direct sequence analysis of exon 1 of the 4 patients with an aberrant MLPA pattern for exon 1 was undertaken to exclude the presence of a SNP interfering with the MLPA probes. No abnormality was identified, indicating that the MLPA results were very likely due to deletions of this region. To delineate all deleted regions, Q-PCR analyses were performed at several points upstream and downstream of the exon(s) involved in the deletions, followed by LR-PCR using the Q-PCR primers mapping just outside the deleted regions. The breakpoints were identified by sequencing the aberrant LR-PCR products. All 4 exon 1 deletions had different breakpoints, all resulting in a complete loss of exon 1 (Figure 2). Three of the 4 deletions did not show a deletion only. In patients 30628 and 21722 also an inversion of 56 and 156 nucleotides, respectively, was present. In patient 21899 an even more complex rearrangement of two deleted regions and an insertion of 7 nucleotides was identified. The exon 1 deletion observed in patient 1264 and the multi-exon deletions in the other families (31457, 29445, 28121 and 14249) did not contain inserted or inverted nucleotides. Only the deletion in patient 28121 was entirely intragenic. The 5' breakpoint was located in intron 8 and the 3' breakpoint in the 3' UTR of exon 23.

The three other deletions started either in the TSC1 upstream region, in intron 1 or in intron 8 and extended into the TSC1 downstream region. None of the 8 breakpoint junctions showed a sequence that could be an obvious trigger for the rearrangements (Table 2). The TSC1 promoter region is located between nucleotide positions 16271 and 15683 upstream of the ATG codon in exon 3 [START_REF] Ali | Identification of a core promoter and a novel isoform of the human TSC1 gene transcript and structural comparison with mouse homolog[END_REF] , indicating that 3 out of the four patients with a complete exon 1 deletion also lack the promoter region. Patient 1264 had a partial deletion of the promoter region. The 155 nucleotides of the 5' end of the promoter region were still present in this patient.

RNA analysis

Because exon 1 is a noncoding exon, it was not clear whether deletion of this exon would be pathogenic, nor whether it would have an effect on the expression of TSC1. Only one patient (21899) was heterozygous for a coding SNP in TSC1 (c.2829C>T in TSC1 exon 22 (rs4962081) (Figure 3A), allowing to assess which alleles of TSC1 were expressed in this patient. To demonstrate equal expression of both alleles of TSC1 in cultured fibroblast cells, control RNA from another individual heterozygous for this SNP, was analysed by RT-PCR. Expression of both alleles could be demonstrated in the control RNA, since the RNA showed a heterozygous pattern for the SNP (Figure 3B). In contrast, RNA from patient 21899 showed only the T nucleotide of SNP rs4962081, indicating monoallelic expression of TSC1 (Figure 3C). We concluded that the deletion of exon 1 in this patient prevented TSC1 expression and deletions affecting this noncoding exon are therefore pathogenic mutations.

Deletion specific PCR analysis

Because the MLPA tests are sensitive to the quality of the DNA [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] (AvdO, unpublished observations), deletion specific PCRs were developed for diagnostic application of mutation analysis within the family, including prenatal testing. In all cases, 3 primers were used: one common primer, one primer located in the deleted region and one primer just outside the deletion (Table 1). The primers were chosen so that the fragment specific for the deletion was shorter than the wildtype fragment (Figure 4). Using the deletion specific PCRs, the unaffected parents of patients 30628, 21722 and 14249 tested negative for the deletion, suggesting that these patients are de novo TSC cases. Patients 1264 and 28121 each had an affected family member. They were available for DNA deletion analysis and tested positive for the respective deletion specific PCRs (Figure 4). The affected sib of patient 1264 showed mild mental retardation, epilepsy, cortical tubers and subependymal nodules (age of diagnosis 15 yrs). The mother of patient 28121 had facial angiofibroma, ungual fibroma, fibrous plaques, hypomelanotic macules and shagreen patches (diagnosed at age 31 yrs.). None of the healthy relatives tested positive for the familial deletions and, thus, nonpenetrance was not encountered in these families. The parents of patient 31457 were tested elsewhere by MLPA analysis. Both parents showed a normal result (data not shown). The parents of individuals 21899 and 29445 were not available for testing.

Clinical details

With one exception [START_REF] Nellist | Large deletion at the TSC1 locus in a family with tuberous sclerosis complex[END_REF] , no clinical information was available on the other previously reported TSC1 deletion cases. The clinical features of the 9 TSC1 deletion index patients are summarised in Supplementary Table 3. All patients had a definite clinical diagnosis of TSC. The clinical findings of these 9 TSC patients with TSC1 deletions were compared to other patients with a TSC1 mutation [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] . Although the number of patients was very small, making comparisons difficult, fewer neurological symptoms and dermatological findings, especially shagreen patches, were found in the TSC1 deletion patient group (Supplementary Table 3; compare last 2 columns).

DISCUSSION

In TSC patients, different types of mutations can be identified in TSC1 or TSC2. Approximately one third are identified in TSC1 [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF][START_REF] Dabora | Mutational analysis in a cohort of 224 tuberous sclerosis patients indicates increased severity of TSC2, compared with TSC1, disease in multiple organs[END_REF] . Most of the pathogenic mutations are nonsense, frameshift or splice site mutations and some missense mutations have been described (http://chromium.liacs.nl/LOVD2/TSC/home.php).

Recently, functional tests have helped to classify missense changes as pathogenic mutations in TSC1 [START_REF] Coevoets | A reliable cell-based assay for testing unclassified TSC2 gene variants[END_REF][START_REF] Nellist | Functional characterisation of the TSC1-TSC2 complex to assess multiple TSC2 variants identified in single families affected by tuberous sclerosis complex[END_REF][START_REF] Nellist | Missense mutations to the TSC1 gene cause tuberous sclerosis complex[END_REF] . So far, only a small number of large rearrangements in TSC1 have been described [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF][START_REF] Longa | TSC1 and TSC2 deletions differ in size, preference for recombinatorial sequences, and location within the gene[END_REF][START_REF] Nellist | Large deletion at the TSC1 locus in a family with tuberous sclerosis complex[END_REF] . In total, 16 large rearrangements, including the deletions described here, have been identified. In our cohort of individuals with a TSC2 mutation (487), an MLPA abnormality was present in 48 cases (9.9%; data not shown). In one case, an intragenic duplication of several exons was identified, whereas the other cases had (multiple) exon deletions. A TSC1 pathogenic mutation was present in 181 individuals (data not shown). Of these mutations, 9 were large rearrangements (5.0%). The percentage large TSC1 rearrangements in the patient group of Kozlowski et al. was lower (0.5% of all TSC1 mutations) compared to our group [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] . In our cohort, all TSC1 and TSC2 rearrangements (n= 57) account for 8.5% of all mutations (n= 668). This is comparable with the data of Kozlowski et al. (6.1% of all TSC mutations).

Of all 16 patients/families with a large TSC1 deletion, 5 patients were sporadic, 5 were familial cases and of the remaining 6 cases, the parents were not analysed by molecular techniques. In our previously described cohort of patients with point mutations in TSC1, a comparable ratio of sporadic to familial cases was observed (22 sporadic to 20 familial) [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] .

Since exon 1 is a noncoding exon, RNA analysis was performed to investigate the effect of deleting this exon on TSC1 mRNA expression. Monoallelic expression was demonstrated, indicating that exon 1 deletions are likely to be null alleles and that there is no alternative promoter present in the TSC1 region. Regulatory elements necessary for basal transcription of TSC1 and a region for optimal promoter function have been defined and were localised to the regions between nucleotide positions 16271 -16003 and 16002 -15683 respectively, upstream of the start codon in exon 3 [START_REF] Ali | Identification of a core promoter and a novel isoform of the human TSC1 gene transcript and structural comparison with mouse homolog[END_REF] . In the UCSC Genome Browser only one promoter region is presented. A CpG island containing region and several transcription binding sites are located in the same region as defined by the functional test and the monoallelic TSC1 expression in our patient. The promoter region was completely deleted in patients 30628, 21722 and 21899. In patient 1264, a partial deletion of the promoter region was identified; only the most 5' 155 nucleotides of the "basal" transcription core were present.

Unfortunately, this patient was not heterozygous for a coding SNP and therefore it was not possible to analyse monoallelic expression of TSC1.

Three of the 4 exon 1 abnormalities were complex events (Figure 1 and Table 2). A combination of a deletion and an inversion was detected in 2 patients, whereas the abnormality in the third patient consisted of 2 deletions separated by 45 nucleotides and an insertion of unknown origin. It is not clear why exon 1 deletions account for almost half of the large TSC1 rearrangements in our cohort. Most deletions were not associated with specific repeat sequences and in only 2 cases a very short repeat sequence (2 or 3 nucleotides) was observed. This is in contrast to the large rearrangements described in TSC2 [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] , where in 70% of the abnormalities very short sequence repeats were present at the junction of the deleted segments.

Although genotype-phenotype comparisons with such a small number of TSC patients should be made with caution, the clinical phenotype of the patients with a TSC1 deletion was slightly less severe overall than that of patients with other TSC1 mutations. In addition, we noted that all patients with a deletion of exon 1 had epilepsy, while this was only observed in 1 of 5 individuals with a deletion affecting other exons of TSC1. We demonstrated that deletions encompassing exon 1 are true null alleles. Therefore, one possible explanation for our observation is that the expression of truncated or mutant TSC1 isoforms may modify the TSC phenotype. Mutant TSC1 isoforms could either have a dominant negative effect by competing with wild-type TSC1 to form inactive TSC1-TSC2 complexes, or have a protective effect by retaining some functionality and maintaining some TSC1-TSC2 activity in the cell [START_REF] Mozaffari | Identification of a region required for TSC1 stability by functional analysis of TSC1 missense mutations found in individuals with tuberous sclerosis complex[END_REF] (MN unpublished observations).

In our cohort, 5.0% of all TSC1 mutations were identified using MLPA, indicating that it is necessary to screen for large TSC1 rearrangements in TSC patients. Although the number of patients identified with a large (complex) deletion is relatively small, it might be that these patients show less severe symptoms compared to patients with point mutations in TSC1. given [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] . On the right side of each deletion the identification number of the patient is given, followed by the indication S (sporadic), F (Familial), ND (not determined) and, if previously published, the references: *Kozlowski et al. [START_REF] Kozlowski | Identification of 54 large deletions/duplications in TSC1 and TSC2 using MLPA, and genotypephenotype correlations[END_REF] , **Longa et al. [START_REF] Longa | TSC1 and TSC2 deletions differ in size, preference for recombinatorial sequences, and location within the gene[END_REF] , ***Nellist et al. [START_REF] Nellist | Large deletion at the TSC1 locus in a family with tuberous sclerosis complex[END_REF] .

Supplementary Table 3. Clinical features of the index patients

The identification numbers for the patients with TSC1 deletions identified in our cohort and of family 1219 of which clinical data were published [START_REF] Nellist | Large deletion at the TSC1 locus in a family with tuberous sclerosis complex[END_REF] are given in the top row. The clinical data from these families were compared to the overall TSC1 mutation positive families previously described [START_REF] Sancak | Mutational analysis of the TSC1 and TSC2 genes in a diagnostic setting: genotype--phenotype correlations and comparison of diagnostic DNA techniques in Tuberous Sclerosis Complex[END_REF] . Gender: M, male and F: female. + Clinical feature present; 0 absent; NE: not examined; ND: no data. If the age of diagnosis is not known, the age of diagnosis is calculated by subtracting the year of birth of the patient of the year the laboratory received the request for mutation analysis.

Figure 1 .

 1 Figure 1. MLPA results. Shown are the graphs after analysis with the Genemarker

Figure 2 .

 2 Figure 2. Overview of the TSC1 deletions identified during this study and described

Figure 3 .

 3 Figure 3. RNA analysis of the coding SNP rs4962081 in exon 22 of TSC1.

Figure 4 .

 4 Figure 4. Agarose gel electrophoresis of the deletion PCR products. The identification

  Figure 1

Table 1

 1 Oligonucleotides used in this study.

	SNP for Q-PCR analysis	Sequence	Ref. seq. NM_000368.3
		5'-3'	5'-3'
	rs36021960	ttcctgtcttgctcggttactca	134.838.750
		ttttcggtgagactggg	134.838.712
		catgacgttcctgcccttaag	134.838.670
	rs2905078	ggctggaagtcggaaatcaa	134.831.005
		tgtcggcaggattg	134.830.983
		cctcggagactccagaaggaa	134.830.948
	rs3011289	cctgccaaagtacagcagtttg	134.821.751
		ctccaagtgcaccctt	134.821.773
		tccagcaacaggtggtacattt	134.821.812
	rs7040593	acccttgacagtggaggacatt	134.816.600
		cgagaactcttcatcgac	134.816.626
		ttgcagttccaagagaggtttga	134.816.670
	rs869116	ctccctccgcaccagatg	134.814.011
		agcctgtgctggtca	134.814.037
		cccagacagagaagggcaaa	134.814.078
	rs12380834	agttccaccaatctgcaaactacttt	134.811.775
		aggaacagactcttccct	134.811.803
		gaatcccttccctggtgtgaa	134.811.846
	rs11796704	ccagccttcttgtttccataaaat	134.807.844
		agtgcttcaggtcctg	134.807.817
		ccagtgaggaagaggacaactga	134.807.775
	rs700797	tggactcagttgccctctgaa	134.805.668
		ctgtgggaagctatgg	134.805.644
		tgggacataaagggtagaagagaaa	134.805.598
	rs13300390	gctctctcacagctcataatgcat	134.803.943
		ctagggcagagatatgtgaa	134.803.914
		caggacccgtacgccattt	134.803.874
	rs5203101	ttaggaggagccaaaggtagactct	134.801.009
		taggctcaggaatggg	134.800.982
		ccttggctaagccacatgct	134.800.946
	rs12337302	ctggttctcgtctgtgcctagtac	134.784.975
		cagctgtcatcctagtct	134.784.950
		cacaggagagaggcgaagga	134.784.903
	rs2809244	tcactggctccttcctaccaa	134.760.050
		ccctgcctcagctg	134.760.075
		gcctaagaactgtggtctggtgtt	134.760.113
	rs2073869	tcctgcagataccctcatgatg	134.753.491
		tcagctgtgacgaggc	134.753.524
		acgccgccgtagtggtt	134.753.560
	rs2231405	tgcttggcatcccacagtt	134.749.724
		ctcaacccccagtgga	134.749.704
		ggatcctgtccttcccatca	134.749.666
	rs2519759	agagggaaaatggcacagtca	134.744.564
		tgcagacagcccc	134.744.586
		aagcggagccagaacttgaa	134.744.621
	rs2072058	tgaaacccaacaggcgacat	134.739.974
		agggcaggcgaaag	134.739.999
		cagaatgacccatggaatcaatta	134.740.042
	rs12555164	cacccaaacacaaagagtgtaacg	134.693.461
		tgggacagagcacg	134.693.486
		caggcagaactttccgtcatg	134.693.523
	rs35602700	tggcaactaaggctgcatga	134.673.378
		ttaagctcgaccctgttt	134.673.399
		caggagtttgctctgcatacga	134.673.439
	rs4962216	tgaattcatttcttctcacattttcc	134.654.970
		cccccatctttgagc	134.655.000
		ttctaagcatgacacagcattggt	134.655.039
	Deletion specific PCR
	c.[-24893_-15354del9540; -24902_ -24847inv56]	agctccttgggaaacaggat	134.810.337
		gggctctgcaaatagctgtc	134.809.122
		aaacagggggcaggaaatag	134.819.205

.

ACKNOWLEDGMENTS

The authors would like to thank the members of the families for providing samples and cooperating with this study. Financial support for MN was provided by the US Department of Defense Congressionally-Directed Medical Research Program (grant no TS060052).